Kirjoittaja Aihe: normaalijakauma pythonilla  (Luettu 3976 kertaa)

seniori

  • Käyttäjä
  • Viestejä: 83
    • Profiili
normaalijakauma pythonilla
« : 27.12.07 - klo:09.50 »
Teen pientä python-ohjelmaa (suuria en osaisikaan), jossa tarvitsen
funktiota normaalijakauman kertymäfunktion arvoille. En kuitenkaan
löytänyt help-ikkunasta tällaista. Varmaan sellainen on olemassa,
mutta mistähän sen löytäisi?
Kiitos vastauksesta, mainitsemasi paketit taitavat olla aika laajoja.
En taida osata käyttää niitä, mutta ongelma ratkesi kun tein
pienen funktion, joka laskee halutun asian kuuden merkitsevän numeron tarkkuudella ja
ongelmaani riittää jo neljän numeron tarkkuus.
Katselen vielä joskus numeric pythonin ohjeita, kunhan olisi aikaa.
« Viimeksi muokattu: 27.12.07 - klo:18.32 kirjoittanut seniori »

SuperOscar

  • Käyttäjä
  • Viestejä: 4063
  • Ocatarinetabellatsumtsum!
    • Profiili
    • Legisign.org
Vs: normaalijakauma pythonilla
« Vastaus #1 : 27.12.07 - klo:13.43 »
En löydä minäkään valmista funktiota Pythonin peruskirjastosta, mutta ehkä NumPy- tai SciPy-paketeissa on?
pöytäkone 1, NUC: openSUSE Leap 15.6, kannettavat 1–3: Debian GNU/Linux 12; pöytäkone 2: openSUSE Tumbleweed; RPi 1: FreeBSD 14-RELEASE; RPi 2: LibreELEC 11

Regel

  • Käyttäjä
  • Viestejä: 1090
  • Lucid
    • Profiili
Vs: normaalijakauma pythonilla
« Vastaus #2 : 28.12.07 - klo:10.26 »
Tuonhan voi itsekin kehitellä:

integraali tiheysfunktiosta, laskee vaikka jollain numeerisella menetelmällä, jollei pythonissa ole valmiiksi integrointia.

edit: scipyssä on numeerinen integrointi

Edit2: Eihän normaalijakauman kertymäfunktiolle ole edes alkeisfunktioiden avulla ilmoitettavaa muotoa.
Edit3: No tässä ois valmis moduuli, jos joku vielä tarttee joskus:

s-alkuiset ovat normitettuja normaalijakaumia, muut eivät

Koodia: [Valitse]
#!/usr/bin/python
# -*- coding: ISO-8859-15 -*-
# Filename: normaldistribution.py
# Developer: Janne Solanpää
# Email: janne.solanpaa[ the funny at-sign ]gmail.com
# Depends on python-scipy and python-math
'''     This program is free software: you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program.  If not, see <http://www.gnu.org/licenses/>.'''


from scipy.integrate import quad
from math import pi, sqrt, e

# Standard normal distribution

#Probability density function
def snormaldist_pdf(x):
return 1/sqrt(2*pi)*e**(-0.5*x**2)

#Cumulative distribution function
def snormaldist_cdf(a,b):
return quad(lambda x: 1/sqrt(2*pi)*e**(-0.5*x**2),a,b)[0]

# The non-standard normal distribution

#Probability density function
# o = expected value
# h = standard deviation
def normaldist_pdf(x,o,h):
return 1/(h*sqrt(2*pi))*e**(-((x-o)**2)/(2*h**2))

#Cumulative distribution function
# o = expected value
# h = standard deviation
def normaldist_cdf(a,b,o,h):
return quad(lambda x: 1/(h*sqrt(2*pi))*e**(-((x-o)**2)/(2*h**2)),a,b)[0]

Edit4: Tiivistetty koodia, ja lisätty lisenssi sekä selitykset.
« Viimeksi muokattu: 28.12.07 - klo:11.37 kirjoittanut Regel »

Heikki Mäntysaari

  • Käyttäjä / tiedottaja
  • Viestejä: 377
    • Profiili
Vs: normaalijakauma pythonilla
« Vastaus #3 : 28.12.07 - klo:11.12 »
Tuonhan voi itsekin kehitellä:
Tiheysfunktiolla ei tosiaan ole alkeisfunktioiden avulla esitettävää integraalifunktiota, joten kertymä on laskettava numeerisesti. Tässä esimerkki laskusta käyttäen numeeriseen integrointiin Simpsonin sääntöä. Alkuparametreja säätämällä voi vaikuttaa nopeuteen ja tarkkuuteen:
Koodia: [Valitse]
import math

# Jakovälien määrä, oltava parillinen!
jakovaleja = 500000
# Mistä aloitetaan numeerinen integrointi, jokin pieni luku
alku=-7

def njak(x):
    return 1/math.sqrt(2*math.pi)*math.exp(-0.5*pow(x,2))


def kertyma(x):
    # Simpsonin säännöllä
    sum=0
    #ensimmäinen termi
    sum+=njak(alku)

    # Jakovälin pituus
    pituus = abs(x-alku)/jakovaleja

    for i in range(1,jakovaleja-1):
        if (i%2==1):
            sum+=4*njak(alku+i*pituus)
        if (i%2==0):
            sum+=2*njak(alku+i*pituus)

        #    print "sum: %f , i: %i , x: %f,  njak: %f" %( sum, i, alku+i*pituus, njak(alku+i*pituus))

    # Viimeinen termi
    sum+=njak(alku+jakovaleja*pituus)

    # Lopuksi
    sum*=(pituus/3.0)

    return sum

print kertyma(0.41)
Suomenkielinen Linux-wiki: Linux.fi - katso myös http://linux.fi/foorumi

Regel

  • Käyttäjä
  • Viestejä: 1090
  • Lucid
    • Profiili
Vs: normaalijakauma pythonilla
« Vastaus #4 : 28.12.07 - klo:11.39 »
Mitenkähän pythonilla pyöristetään?

Meinaan, että tulee välillä hieman päälle 1 kokonaistodennäköisyydeksi tuolla ohjelmallani. Voisi vähän viimeisiä desimaaleja pätkiä.

seniori

  • Käyttäjä
  • Viestejä: 83
    • Profiili
Vs: normaalijakauma pythonilla (ratkaistu)
« Vastaus #5 : 28.12.07 - klo:20.36 »
Kiitokset Heikille ja Regelille. Heikin ohjelman tarkkuus on yhtä vaatimaton kuin oman
ohjelmanikin mutta regelin scipyyn perustuva ohjelma on tarkka
ja vieläpä tosi nopea. "Tarkistin" sen laskemat arvot openofficen calcissa ja
hyvältä näyttää. Tätä scipyä tulen käyttämään, KIITOKSET.
Calcissa:
NORMDIST(1,0,1,1) jossa
1. luku on sat. muuttujan lukuarvo,
2. luku on sen jakauman odotusarvo(keskiarvo)
3. luku on keskihajonta (STDEV) ja
4. luku kertoo lasketaanko tiheysfunktion arvo (C=0) vai
lasketaanko kertymäfunktion arvo (C=1).Tulos 0,34144...
PS. Olisi mukava tietää perustuuko openofficen menetelmä esim.
rombergin metodiin vai johonkin muuhun, simpsonin sääntöön
se ei varmaankaan perustu.
Onkohan tällaisesta missään saatavissa tietoa?

Regel

  • Käyttäjä
  • Viestejä: 1090
  • Lucid
    • Profiili
Vs: normaalijakauma pythonilla
« Vastaus #6 : 28.12.07 - klo:22.44 »
Mukava olla avuksi  :)

Tulipa vielä itsekin löydettyä tuo scipy  :D Hieno moduuli, pakko myöntää.