Doç. Dr. Özgür Baştürk
Ankara Üniversitesi, Astronomi ve Uzay Bilimleri Bölümü
obasturk at ankara.edu.tr
http://ozgur.astrotux.org
Olasılık Dağılım Fonksiyonlarının Özellikleri¶
Rastgele Değişken (ing. random variable): Olası değerleri rastgele bir olayın sayısal sonuçları olan bir değişkendir. İki tür rastgele değişken vardır, süreksiz (kesikli) (ing. discrete) ve sürekli (ing. continuous) değişkenler.
Süreksiz bir rastgele değişken, sadece sayılabilir değerler alabilen ve böylece nicelleştirilebilen bir değişkendir. Örneğin, adil bir zarın her atılıştaki değeri rastgele bir $x$ değişkeni olarak tanımlanabilir. Bu durumda $x$, $[1,2,3,4,5,6]$ değerlerini alabilir. Bu nedenle de $x$, süreksiz rastgele bir değişkendir.
Süreksiz rastgele bir değişkenin olasılık dağılımı, olası değerlerinin her biri ile ilişkili olasılıkların bir listesidir. Buna olasılık fonksiyonu ya da daha doğru bir ifadeyle Olasılık Kütle Fonksiyonu (ing. Probability Mass Function, PMF) denir.
Rastgele bir $X$ değişkeninin k farklı değer alabileceğini varsayalım. $X = x_i$ 'ye karşılık gelen $P (X = x_i) = p_i$ olasılıkları aşağıdaki koşulları sağlar:
- $0 \lt p_i \lt 1$, $i = 1, 2, 3, ...$
- $p_1 + p_2 + ... + p_k = 1$
Sürekli rastgele değişken, sonsuz sayıda olası değeri alabilen değişkendir. Örneğin, bir sınıftaki öğrencilerin boyları rastgele bir $x$ değişkeni olarak tanımlandığında sürekli bir rastgele değişkendir Sürekli rastgele değişken bir değerler aralığı içinde tanımlandığından, o aralıkta olma olasılığı, eğrinin (veya integralin) altındaki alanla temsil edilir.
Olasılık Yoğunluk Fonksiyonu (ing. Probability Density Function, PDF olarak bilinen sürekli rasgele değişkenin olasılık dağılımı, sürekli değerleri alan fonksiyonlardır. Herhangi bir tek değeri gözlemleme olasılığı, 0'a eşittir, çünkü rastgele değişken tarafından alınabilecek değerlerin sayısı sonsuzdur ($p(x_i) = 0$. Olasılık yoğunluk fonksiyonu için $P (a \lt x \lt b) = p$ olasılıkları aşağıdaki koşulları sağlar:
Olasılık her zaman pozitiftir ($x$'in bütün değerleri için $p \gt 0$).
$p(x)$'in altında kalan 1'dir.
Ayrıca, tüm rastgele değişkenlerin (sürekli ya da değil) birikimli (kümülatif) dağılım karşılğı da vardır. Birikimli dağılım fonksiyonu, her $x$ değeri için rastgele değişken $X$'in $x$'den küçük veya ona eşit olma olasılığını veren bir fonksiyondur ($P(X \le x) = p(x))$. Kesikli rastgele bir değişken için olasılık dağılımları toplanarak birikimli (kümülatif) dağılım fonksiyonu bulunur.
Python'da Dağılımlar¶
Python'da tüm olasılık dağılım fonksiyonlarına ilişkin fonksiyonların yer aldığı modüller bulunmaktadır. Bu modüllerden en sık kullanılanları scipy.stats ve np.random paketleridir. Bu modüllerdeki fonksiyonlarla ilgilenilen dağılımın yapısına uygun örnekler ya da dağılım fonksiyonları oluşturularak incelenebilir. Çizim için matplotlib.pyplot fonksiyonları hist, bar ve plot en sık kullanılanlardır. Ayrıca seaborn modülü fonksiyonları da kolay kullanılabilmeleri nedeniyle önemli bir alternatif oluşturmaktadır.
scipy.stats paketindeki neredeyse tüm dağılım fonksiyonlarında bulunan bazı argüman ve parametreler, Normal Dağılım örneğinde (stats.norm) aşağıda listelenmeye çalışılmıştır.
rvs(loc=0, scale=1, size=1, random_state=None): Söz konusu dağılımdan bir örnek oluşturmak için kullanılır.
pdf(x, loc=0, scale=1): Olasılık yoğunluk fonksiyonu oluşturmak için kullanılır.
logpdf(x, loc=0, scale=1): Olasılk yoğunluk fonksiyonun logaritması (doğal logaritma, ln)'nı verir.
cdf(x, loc=0, scale=1): İlgilenilen dağılım türü temelinde bir birikimli (kümülatif) yoğunluk fonksiyonu oluşturmak için kullanılır.
logcdf(x, loc=0, scale=1): Birikimli yoğunluk fonksiyonunun logaritmasını (ln) verir.
sf(x, loc=0, scale=1): Kümülatif fonksiyon $P(X \le x)$ olasılığını verirken $sf$ (ing. survival function) $P(X \gt x) = 1 - (P(X \le x)$ olasılığını verir.
logsf(x, loc=0, scale=1): sf fonksiyonunun logaritmasıdır.
ppf(q, loc=0, scale=1): (ing. percent point function), birikimli yoğunluk fonksiyonunda elde edilen olasılığın karşılık geldiği $X = x$ değeridir. Bu anlamda birikimli yoğunluk fonksiyonunun tersidir. Birikimli yoğunluk fonksiyonu $X = x$ değeri argüman olarak sağlandığında $P(X \ge x)$ olasılığını verirken, $ppf$ verilen $P(X \ge x)$ olasılık değeri için $X = x$ değerini verir.
isf(q, loc=0, scale=1): (ing. inverse survival function), ppf'e benzer şekilde $1 - cdf$ olasılığının karşılığı olan $X = x$ değerini verir.
moment(n, loc=0, scale=1): Dağılımın n. momentini verir.
stats(loc=0, scale=1, moments=’mv’): Dağılıma ilişkin momentleri verir. İstenen momentler
momentsparametresine baş harfiyle sağlanır. (Ortalama: Mean(‘m’), Varyans: variance(‘v’), Çarpıklık: skew(‘s’), ve/veya kurtosis(‘k’)).entropy(loc=0, scale=1): Dağılımın entropisini verir.
fit(data, loc=0, scale=1):
dataparametresi üzerinden verilen bir veriye ilgili dağılımı uyumlar (fit eder) ve bu uyumalamanın dağılımın türüne göre ilgili parametrelerini döndürür.expect(func, args=(), loc=0, scale=1, lb=None, ub=None, conditional=False, **kwds): Tek değişkenli bir fonksiyonun bir dağılıma göre beklenen değerini döndüren fonksiyondur.
median(loc=0, scale=1): Dağılımın ortancasını (medyan değerini) döndürür.
mean(loc=0, scale=1): Dağılımın ortalama değerini döndürür.
var(loc=0, scale=1): Dağılımın varyansını ($\sigma^2$) döndürür.
std(loc=0, scale=1): Dağılımın standart sapmasını ($\sigma$) döndürür.
interval(alpha, loc=0, scale=1): Tüm dağılımın $\alpha$ yüzdelik bölümüne denk gelen aralığının uç noktalarını verir.
Tekdüze Dağılımlar¶
Bir değişkenin tekdüze dağılımı onun bir aralık içerisindeki değerlerin tamamının aynı olasılıkla alabileceği sürekli bir dağılımdır.
numpy.random fonksiyonlarından uniform ile tekdüze bir dağılımdan giderek daha fazla sayıda örnek seçilerek tekdüze bir dağılımın nasıl oluştuğunu örnekleyen bir dizi kod aşağıda verilmiştir. Ayrıca tekdüze dağılımın her aşamada bir normal fonksiyona göre ne kadar daha çarpık (ing. skewness) ve basık olduğuna ilişkin değerler (dağılımın $skew$ ve $kurtosis$ parametreleriyle) de verilmiştir. Tekdüze bir dağılımın simetrik olması beklendiği için örnek sayısı arttırıldıkça çarpıklığının 0'a, ortancaya (medyan), basıklığı gösteren Kurtosis değerinin ise -1.2'ye yaklaşması beklenmelidir.
import numpy as np
from matplotlib import pyplot as plt
from scipy import stats as st
%matplotlib inline
x = np.random.uniform(-5,5,10)
print(x.mean(), x.std(), np.median(x))
plt.hist(x, bins=20)
plt.show()
x = np.random.uniform(-5,5,100)
print(x.mean(), x.std(), np.median(x))
plt.hist(x, bins=20)
plt.axvline(x = x.mean(), ls='--', c='r')
plt.axvline(x = np.median(x), ls='--', c='g')
plt.show()
print('Normal bir dagilima gore basiklik: ', st.kurtosis(x))
print('Normal bir dagilima gore carpiklik: ',st.skew(x))
x = np.random.uniform(-5,5,1000)
print(x.mean(), x.std(), np.median(x))
plt.hist(x, bins=20)
plt.axvline(x = x.mean(), ls='--', c='r')
plt.axvline(x = np.median(x), ls='--', c='g')
plt.show()
print('Normal bir dagilima gore basiklik: ', st.kurtosis(x))
print('Normal bir dagilima gore carpiklik: ',st.skew(x))
x = np.random.uniform(-5,5,100000)
print(x.mean(), x.std(), np.median(x))
plt.hist(x, bins=20)
plt.axvline(x = x.mean(), ls='--', c='r')
plt.axvline(x = np.median(x), ls='--', c='g')
plt.show()
print('Normal bir dagilima gore basiklik: ', st.kurtosis(x))
print('Normal bir dagilima gore carpiklik: ',st.skew(x))
N = 100000
x = np.random.uniform(-5,5,N)
print(x.mean(), x.std(), np.median(x))
plt.xlabel('x')
plt.ylabel('PDF')
plt.hist(x, bins=20, density=True)
plt.axvline(x = x.mean(), ls='--', c='r')
plt.axvline(x = np.median(x), ls='--', c='g')
plt.show()
N = 100000
x = np.random.uniform(-5,5,N)
print(x.mean(), x.std(), np.median(x))
plt.xlabel('x')
plt.ylabel('PDF')
plt.hist(x, bins=20, density=True)
plt.axvline(x = x.mean(), ls='--', c='r')
plt.axvline(x = np.median(x), ls='--', c='g')
plt.show()
Normal Dağılım ve Gauss Fonksiyonu¶
Gauss dağılımı olarak da bilinen Normal Dağılım, deneysel bütün bilimlerde yeri bulunan sürekli bir olasılık dağılımıdır. Normal dağılım, ortalama, $\mu$ ve standart sapma, $\sigma$) ile tanımlanan çan şeklindeki yoğunluk eğrisine sahiptir. Yoğunluk eğrisi (ing. Probability Density Function, PDF), merkezi ortalamaya yakın, yayılımı standart sapması ile belirlenen; ortalamaya yakın verilerin, ortalamadan çok uzaktaki verilere göre daha sık olduğunu gösteren simetrik bir fonksiyondur. Belirli bir $x$ noktasında ortalama $\mu$ ve standart sapma $\sigma$ ile normal yoğunluk eğrisinin olasılık dağılım fonksiyonu şu şekilde verilir:
$$ p(x)~dx = \frac{1}{\sigma\sqrt{2\pi}} e^{\frac{-1}{2} (\frac{x - \mu}{\sigma})^2}~dx $$def gauss(x,mu,sigma):
a = 1/(np.sqrt(2*np.pi)*sigma)
b = np.exp(-0.5*((x - mu)/sigma)**2)
return a*b
x = np.random.randn(100000)
print(x.mean(), x.std(), np.median(x))
y = np.linspace(-5,5,100000)
fy = gauss(y,x.mean(), x.std())
plt.hist(x, bins=50, density=True)
plt.plot(y,fy,'r-')
plt.axvline(x = x.mean(), ls='--', c='r')
plt.axvline(x = np.median(x), ls='--', c='g')
plt.show()
print('Normal bir dagilima gore basiklik: ', st.kurtosis(x))
print('Normal bir dagilima gore carpiklik: ',st.skew(x))
Örnek: Örneğin RR Lyrae yıldızlarının zonklama dönemlerini inceliyor olalım. Evrendeki bütün RR Lyrae yıldızlarını gözleyemeyeceğimiz için RR Lyrae değişenleri popülasyonundan ancak örnekler alabiliriz.
Bugüne kadar yapılan gözlemler bize RR Lyrae yıldızlarının zonklama dönemlerinin 4 saat ile 24 saat arasında değiştiğini göstermektedir. Ortalama zonklama döneminin 14 saat civarında olduğunu ve standart sapmasının da 5 saat olduğunu varsayalım.
Port = 14. # ortalama zonklama donemi
Pstd = 5 # zonkalama doneminin standart sapmasi
Şimdi 500 ayrı küresel kümede 20'şer RR Lyrae yıldızının gözlendiğini varsayalım ve bu örneklerin ortalamasına ve standart sapmasına bakalım. Bunun için her biri 20 yıldız içeren 500 grup için zonklama döneminin normal dağıldığını varsayarak ve bunun için np.random.normal fonksiyonunu kullanarak rastgele zonklama dönemleri verisi oluşturup, bu grupların ortalaması ve standart sapmasına bakalım.
gozlemler = np.zeros((500,20))
for i in range(gozlemler.shape[0]):
for j in range(gozlemler.shape[1]):
gozlemler[i,j] = Port + np.random.normal()*Pstd
print('Tum orneklerin ortalamasi:', round(np.mean(gozlemler),4))
print('Standart sapmasi:', round(np.var(gozlemler)**0.5,4))
Şimdi her bir küresel kümenin ortalamalarının nasıl dağıldığına bir histogramla bakalım.
# Her bir kuresel kumedeki 20ser RR Lyrae
# yildizinin zonklama donemi ortalalamalarinin dagilimi
from matplotlib import pyplot as plt
fig, ax = plt.subplots(figsize=(9,6))
plt.hist(np.mean(gozlemler,axis=1),
bins=20, label='Ortalama Donem')
ax.set_xlabel("Zonklama Donemi (saat)",fontsize=12)
ax.set_ylabel("Frekans",fontsize=12)
# Aslinda tum RR Lyrae'lar icin ortalama
# parlakligi 14 gun varsayiyoruz
plt.axvline(x=14, ls='--', c='r')
plt.legend()
plt.tight_layout()
Görüldüğü üzere her bir küresel kümedeki RR Lyrae yıldızlarının zonkalama dönemlerinin ortalaması da normal dağılmış durumdadır. Herhangi bir küresel kümeyi (örneğin 42.'yi) alsak ve tüm kümelerde RR Lyrae yıldızlarının dağılımıyla karşılaştırsak:
fig, ax = plt.subplots(figsize=(9,6))
# 500 kumede RR Lyrae yildizlarinin
# ortalama zonklama donemlerinin dagilimi
plt.hist(np.mean(gozlemler,axis=1),
bins=20, alpha=0.2, label='Tum Kumeler')
# Rastgele sectigimiz 42. kume
kume42ort = np.mean(gozlemler[41])
kume42std = np.var(gozlemler[41])**0.5
# Secilen kumenin ortalama zonklama doneminin hatasi (standart hata)
std_hata = kume42std/(gozlemler[41].shape[0])**0.5
# Sadece bu kumeden yola cikarak 500 kuresel kume icin
# nasil bir dagilim elde edecegimize bakalim
cikarim = np.zeros((500,20))
for i in range(cikarim.shape[0]):
for j in range(cikarim.shape[1]):
cikarim[i,j] = kume42ort + np.random.normal()*kume42std
# Simdi bu kumeden yola cikarak elde edecegimiz dagilima bakalim
plt.hist(np.mean(cikarim, axis=1), bins=20, alpha=0.5,
label='42. kumeden cikarim')
plt.xlabel("Ortalama Zonklama Donemi (saat)",fontsize=12)
plt.ylabel("Frekans",fontsize=12)
plt.legend()
plt.show()
Sadece 20 RR Lyrae yıldızının ortalama parlaklıkları ve standart hatasını kullanarak 500 küme için oluşturduğumuz dağılım her ne kadar normal bir dağılım olsa da ortalaması ve standart sapması orjinal gözlemlerimizdeki dağılımdan bir miktar farklı oldu. Yine de sadece 20 yıldıza ilişkin istatistikten yola çıkarak oluşturduğumuz bu dağılım, ortalama da normal dağıldığı için orjinal verimizdeki dağılımdan çok farklı değil.
Kümelerin standart sapmalarının dağılımına baktığımızda onun da normal dağıldığını görebiliriz.
# standart sapmalarin dagilimi
gozlemlerstd = np.var(gozlemler, axis=1)**0.5
# standart sapmanin histogrami
fig, ax = plt.subplots(figsize=(9,6))
plt.hist(gozlemlerstd, bins=20, label='Standart sapma')
plt.xlabel("Standart Sapma (saat)",fontsize=14)
plt.ylabel("Frekans",fontsize=14)
# verimizi olustururken std'yi 10 almistik
plt.axvline(x=Pstd,ls="--", c="r")
plt.legend()
plt.tight_layout()
plt.show()
Standart Normal Dağılım¶
Standart normal dağılımda ($\mu = 0$, $\sigma = 1$) herhangi bir $x$ değişkeninin ortalamadan z = $\pm \Delta x$ kadar sapması olasılığı scipy.stats.norm paketi kümülatif dağılım fonksiyonu (ing. cumulative distribution function) cdf kullanılarak elde edilebilir. Matematiksel olarak,
import scipy.stats as st
z = 1.645
pz = st.norm.cdf(z)
print("Standart normal dagilimda x = +{:.2f} sigma sapma olasiligi : {:.2f}". \
format(z,pz))
Fonksiyon $x$'in $-\infty$'dan $\Delta x = +1.65\sigma$ 'ya kadarki aralıkta olma olasılığını verir (standart normal dağılımın bu değere kadar sol tarafının altında kalan alan). Eğer x'in $\Delta x = \pm1.65\sigma$ aralığında olma olasılığı isteniyorsa bu kez bu değerin
z1 = -1.645
z2 = 1.645
pz_z1_z2 = st.norm.cdf(z2) - st.norm.cdf(z1)
print("Standart normal dagilimda x'in = +/-{:.2f} sigma dahilinde olma olasiligi : {:.2f}". \
format(z2,pz_z1_z2))
Dolayısıyla normal dağılımı temsil eden eğrinin en yüksek eğime sahip tanjantlarının eğriyi ve x eksenini kestiği noktalar; sırasıyla 1 ve 2 $\sigma$ verilerek öğrenilebilir.
z1 = -1.0
z2 = 1.0
pz_z1_z2 = st.norm.cdf(z2) - st.norm.cdf(z1)
print("Standart normal dagilimda x'in = +/-{:.2f} sigma dahilinde olma olasiligi : {:.2f}". \
format(z2,pz_z1_z2))
z1 = -2.0
z2 = 2.0
pz_z1_z2 = st.norm.cdf(z2) - st.norm.cdf(z1)
print("Standart normal dagilimda x'in = +/-{:.2f} sigma dahilinde olma olasiligi : {:.2f}". \
format(z2,pz_z1_z2))
Tersten giderek verilen bir olasılığın merkezden kaç standart sapma uzaklığa denk geldiği (yani $z$ değeri) de bulunabilir. Bunun için ppf (percent point function) fonkiyonu kullanılır. Fonksiyonun yine $+\Delta x$'e kadar olasılık verdiği (eğrinin sol tarafının altında kalan) akılda tutulmalıdır.
pz = 0.95
z = st.norm.ppf(pz)
print("Standart normal dagilimda {:.2f} olasiligin karsilik geldigi sapma degeri : {:.2f}". \
format(pz,z))
İki taraflı olasılık isteniyorsa, $px$ argümanı buna uygun olarak verilmelidir.
pz = 0.05
z = st.norm.ppf(pz)
print("Standart normal dagilimda +\-{:.2f} olasiligin karsilik geldigi sapma degeri : +/-{:.2f}". \
format(pz,z))
Aynı amaçla daha kullanışlı olan interval fonksiyonu da kullanılabilir. norm.interval fonksiyonu argüman verilen bir olasılığın $x$'in hangi z = $\pm \Delta x$ değerleri arasında olma olasılığı olduğunu verir. Bir başka deyişle, verilen olasılığa karşılık gelen alanın uç noktaları olan $z$ değerlerini verir.
pz = 0.90
z = st.norm.interval(pz)
print(z)
print("Standart normal dagilimda +\-{:.2f} olasiligin karsilik geldigi sapma degeri : +/-{:.2f}". \
format(pz,z[1]))
import numpy as np
from matplotlib import pyplot as plt
from matplotlib.patches import Polygon
fig, ax = plt.subplots()
x = np.linspace(st.norm.ppf(0.001), st.norm.ppf(0.999), 1000)
ax.plot(x, st.norm.pdf(x), 'r-', lw=2, alpha=0.6, label='norm pdf')
xi = st.norm.ppf(0.001)
xf = st.norm.ppf(0.05)
ix = np.linspace(xi,xf)
iy = st.norm.pdf(ix)
verts = [(xi, 0), *zip(ix, iy), (xf, 0)]
poly = Polygon(verts, facecolor='b')
ax.add_patch(poly)
ax.set_xlabel("z")
ax.set_ylabel("pdf")
plt.arrow(-2.75, 0.10, 0.5, -0.05, capstyle="round",
head_width=0.02, head_length=0.05)
plt.text(-2.6,0.115,"p(z)")
plt.show()
Bu alanın değerini birikmli (kümülatif) dağılım fonksiyonu cdf üzerinden görebiliriz. Bir başka deyişle, kümülatif dağılım fonksiyonu, normal dağılım fonksiyonunun altında kalan alanın z değeri ile değişimini verir. $z \rightarrow \infty$ iken $cdf \rightarrow 1.0$ 'dır.
from matplotlib.patches import Polygon
fig, ax = plt.subplots()
x = np.linspace(st.norm.ppf(0.001), st.norm.ppf(0.999), 1000)
plt.plot(x, st.norm.cdf(x), 'r-', lw=2, alpha=0.6, label='norm pdf')
plt.axvline(x=st.norm.ppf(0.05))
plt.axhline(y=st.norm.cdf(st.norm.ppf(0.05)))
x_degeri = "z = {:.3f}".format(st.norm.ppf(0.05))
y_degeri = "p(z) = {:.3f}".format(st.norm.cdf(st.norm.ppf(0.05)))
plt.text(1, 0.1, y_degeri)
plt.text(-2.8,0.15, x_degeri)
plt.xlabel('z')
plt.ylabel('cdf')
plt.show()
st.norm fonksiyonları sadece standart normal dağılım için değil, herhangi bir normal dağılım için de kullanılabilir. Örneğin RR Lyrae yıldızlarının zonklama dönemleri dağılımından hareketle bir RR Lyrae yıldızının zonklama döneminin 17 günden fazla olma olasılığını bulmak isteyebiliriz. RR Lyrae yıldızları için ortalama dönemin $\mu = 14$ saat, standart sapmasının ise $\sigma = 5$ saat olduğunu belirlemiş ve zonklama döneminin normal dağıldığını varsaymıştık.
Port = 14
Pstd = 5
x = 17
print("Bir RR Lyrae yildizi icin zonklama doneminin {:d} saatten uzun olma olasiligi: {:.2f}".\
format(x, 1 - st.norm.cdf(x, loc=Port, scale=Pstd)))
x = 3
print("Bir RR Lyrae yildizi icin zonklama doneminin {:d} saatten kisa olma olasiligi: {:.2f}".\
format(x, st.norm.cdf(x, loc=Port, scale=Pstd)))
Yukarıdaki ilk örnekte ($P_{zonk} > 17$) görüldüğü gibi cdf normal dağılım eğrisinin sol tarafında kalan alanı verdiği için belirli bir değerden büyük alanı (olasılığı) hesaplamak için $1 - cdf$ kullanmak gerekebilir. Bunun yerine sf fonksiyonu (survival function) da kullanılabilir.
x = 17
print("Bir RR Lyrae yildizi icin zonklama doneminin {:d} saatten uzun olma olasiligi: {:.2f}".\
format(x, st.norm.sf(x, loc=Port, scale=Pstd)))
Aynı şekilde olasılıklardan hareketle zonkalama dönemleri de bulunabilir.
px = 0.95
aralik = st.norm.interval(px,loc=Port, scale=Pstd)
print("RR Lyrae yildizlarinin \%{:.1f}'i {:.2f} ile {:.2f} saatleri arasindadir.".\
format(px*100, aralik[0], aralik[1]))
Log-Normal Dağılım¶
Bir değişimin logaritmasının normal dağılım göstermesi durumunda olasılık yoğunluk fonksiyonu log-normal dağılım fonksiyonu ile temsil edilir. Bir başka deyişle, bir $x$ değişkeni normal dağılıma uygun değerler alıyorsa, $e^x$ log-normal dağılıma uygun değerler alır. log-normal dağılım aşağıdaki şekilde fiade edilir.
$$ p(x)dx = \frac{1}{x~\sigma \sqrt{2\pi}}~exp(-\frac{(ln(x) - \mu)^2}{2\sigma^2})~dx$$Birden fazla pozitif değişkene, o değişkenlerin çarpımıyla bağlı bir değişkenin log-normal dağılım göstermesi beklenir. Eğer $y = exp(x)$ log-normal dağılıyorsa, bu dağılımın ortalama değeri $\bar{y} = exp(\mu + \sigma^2 / 2)$, medianı $exp(\mu)$, modu ise $exp(\mu - \sigma^2)$ ile verilir.
Örnek olarak $\mu = 2$, $\sigma = 1$ parametreleriyle bir log-normal dağılım varsayalım.
import matplotlib.pyplot as plt
import numpy as np
%matplotlib inline
Port = 2.0
Pstd = 1.0
# Ortalamasi Port, standart sapmasi Pstd
# olan log-normal bir dagilimdan alinmis
# 1000 rastgele ornekle bir histogram olusturalim
s = np.random.lognormal(Port, Pstd, 1000)
count, bins, ignored = plt.hist(s, 100, density=True, align='mid')
# histogramin sinirlari dahilindeki 1000 x
# degerini log-normal dagilimin olasilik yogunluk
# fonksiyonunda (pdf) yerine koyalim
# ve histogramin uzerine cizdirelim
x = np.linspace(min(bins), max(bins), 1000)
pdf = 1. / (x * Pstd * np.sqrt(2 * np.pi)) *\
(np.exp(-(np.log(x) - Port)**2 /(2 * Pstd**2)))
plt.plot(x,pdf,'r-')
plt.ylabel('PDF')
plt.axis('tight')
Olasılık yoğunluk fonksiyonunu, dağılımın tanımında (ve kodun içinde $pdf$ değişkeni ile) verilen ifadeyle sağlayabileceğimiz gibi $scipy.stats$ fonksiyonlarından $lognorm$ 'un $pdf$ fonksiyonunu da kullanabiliriz.
import matplotlib.pyplot as plt
import numpy as np
from scipy import stats as st
%matplotlib inline
Port = 2.0
Pstd = 1.0
# Ortalamasi Port, standart sapmasi Pstd
# olan log-normal bir dagilimdan alinmis
# 1000 rastgele ornekle bir histogram olusturalim
s = np.random.lognormal(Port, Pstd, 1000)
count, bins, ignored = plt.hist(s, 100, density=True, align='mid')
# histogramin sinirlari dahilindeki 1000 x
# degerini log-normal dagilimin olasilik yogunluk
# fonksiyonunda (pdf) yerine koyalim
# ve histogramin uzerine cizdirelim
# Bu kez scipy.stats.lognorm kullanalim
x = np.linspace(min(bins), max(bins), 1000)
pdf = st.lognorm.pdf(x,s=Pstd,scale=np.exp(Port))
plt.plot(x,pdf,'r-')
plt.ylabel('PDF')
plt.axis('tight')
Görüldüğü üzere st.lognorm.pdf diğer olasılık yoğunluk fonksiyonlarına göre bir miktar farklı çalışmaktadır. Bir log-normal dağılımın iki parametresi ($\mu$ ve $\sigma$) olmakla birlikte scipy.stats paketindeki diğer dağılım fonksiyonlarında loc dağılımın x-ekseninde kayma, scale parametresi ise ölçeğini verirken; scipy.stats.lognormal fonksiyonlarında $\sigma$, şekil parametresi; $e^{\mu}$, ölçek parametresi olarak verilir. Zira scipy.stats.lognormal fonksiyonları için olasılık yoğunluk fonksiyonu,
olarak tanımlanır. x'in logaritması değil kendisi $\mu$ ve $\sigma$ parametreleri ile normal dağılsaydı;
$$ y = exp(\mu + \sigma ~x) $$$$ y = exp(\mu) ~ exp(x)^\sigma $$$$ y = exp(\mu) ~ exp(x-x_0)^{\sigma} $$olurdu. Bu ifadeyi
$$ Y = m ~ (X - X_0)^n $$şeklinde yazarsak, $m = exp(\mu) = e^{\mu}$'nün ölçek parametresi, $n = \sigma$'nın şekil parametresi $(x - x_0) = e^{X - X_0}$ da $loc = X_0 = ln(x_0)$'dır.
Şimdi ortalamayı sabit tutup standart sapmayı değiştirdiğimiz vakit lognormal dağılımın normal dağılıma nasıl yaklaştığını inceleyelim.
mu = 1.0
sigma = [0.25,0.5,0.75,1.0]
x = np.linspace(0.01, 5.0, 1000)
for sig in sigma:
pdf = (np.exp(-(np.log(x) - mu)**2 / (2 * sig**2))/\
(x * sig * np.sqrt(2 * np.pi)))
leg = "$\sigma$= {:.2f}".format(sig)
plt.plot(x,pdf, label=leg)
plt.xlabel('x')
plt.ylabel('PDF')
plt.grid(True)
plt.legend(loc="upper left")
plt.show()
Aynı grafiği $scipy.stats.lognorm$ fonksiyonlarından $pdf$ ile elde etmeye çalışalım.
mu = 1.0
sigma = [0.25,0.5,0.75,1.0]
x = np.linspace(0.01, 5.0, 1000)
for sig in sigma:
pdf = st.lognorm.pdf(x,s=sig,scale=np.exp(mu))
leg = "$\sigma$= {:.2f}".format(sig)
plt.plot(x,pdf, label=leg)
plt.xlabel('x')
plt.ylabel('PDF')
plt.grid(True)
plt.legend(loc="upper left")
plt.show()
Weibull Dağılımı¶
log-normal dağılımlar gibi asimetrik dağılımlar ailesinin bir üyesi de Weibull dağılımıdır. Başka parametrizasyonları da bulunmakla birlikte olasılık yoğunluk fonksiyonunun (PDF) standart parametriasyonu aşağıdaki gibidir.
$$ f(x, \lambda, k) = \frac{k}{\lambda} (\frac{x}{\lambda})^{k-1}~e^{-(x/\lambda)^{k}} $$$k \gt 0$ olmak üzere şekil parametresi, $\lambda \gt 0$ ölçek parametresidir. X bağımsız parametresi herhangi bir niceliğin azaldığı zamanı (başarısızlığa uğrayana kadarki zamanı) gösteriyorsa, Weibull dağılımı onun zamana kuvvet yasasıyla bağlı olan azalma hızını (ing. failure rate: başarısızlık hızı) veren bir dağılımdır. Bu durumda $k-1$ bu kuvvete karşılık gelir.
Örnek olarak $k \lt 1$ durumu, "başarısızlık hızının" azaldığını gösterir. Örneğin bebeklerdeki "ani bebek ölümü" olarak adlandırılan ölüm olasılığı önemli ölçüde yüksekse bu bebeklerin popülasyondan ayrılmasıyla da hızla azalır. Aynı şekilde herhangi bir şekilde defolu ürünlerin erken bozularak popülasyondan ayrılmaları ürünün başarısız bulunma hızını giderek düşürür. $k = 1$ durumu başarısızlık hızının sabit olduğu duruma karşılık gelir (üstel dağılım). Bu durumda örnekteki bebek ölüm oranlarının ve ürünlerin bozulmalarının dışarıdan bir faktörden etkilendiği düşünülür. $k > 1$ ise başarısızlık oranı zamanla artıyor demektir. Örnekte ürünlerin zamanla yaşlandığı ve çalışmaz hale gelebildiği bir senaryo bu duruma uygundur.
Weibull dağılımı, $k = 1$ durumu üstel dağılıma (ing. exponential distribution), $k = 2$, $\lambda = \sqrt{2} \sigma$ ise Rayleigh dağılımına karşılık gelir.
numpy.random.weibull fonksiyonu bir Weibull dağılımına uyan rastgele bir örnek dizisi oluşturmak için kullanılır.
k = 8. # basarisizlik hizinin zamanla arttigi durum
lamda = 1. # olcek parametresi
wb_ornek = np.random.weibull(k, 1000)
def weib(x,l,k):
return (k / l) * (x / l)**(k - 1) * np.exp(-(x / l)**k)
sayi, bins, _ = plt.hist(np.random.weibull(k,1000))
x = np.arange(1,50.)/25.
olcek = sayi.max()/weib(x, lamda, k).max()
plt.plot(x, weib(x,lamda,k)*olcek)
plt.show()
Tıpkı diğer dağılım fonksiyonlarında olduğu gibi numpy.random fonksiyonlarıyla dağılıma uyan bir örnek oluşturabilir, olasılık yoğunluk fonksiyonu, herhangi bir değerin herhangi bir aralıka olma olasılığı, birikimli dağılım fonksiyonu gibi özellikler ise scipy.stats.weibull_min fonksiyonu ile elde edilebilir ve incelenebilir; dağılıma uygun bir yapı gösterdiği düşünülen örneklere söz konusu dağılım fonksiyonu uyumlanabilir (fit edilebilir). Şimdi $k < 1$ durumu için bir örneğe bakarak "başarısızlık hızının" azalışını görelim.
k = 0.85 # basarisizlik hizinin zamanla azaldigi durum
l = 1. # olcek parametresi
ortalama, varyans, carpiklik, kurtosis = st.weibull_min.stats(k, moments='mvsk')
print("Dagilimin ortalama degeri: {:.2f}, standart sapmasi: {:.2f}, carpikligi: {:.2f} ve basikligi: {:.2f}".\
format(ortalama, varyans, carpiklik, kurtosis))
# Olasilik yogunluk fonksiyonu
# %1 ile %99 oalsilik arasini alalim. Tum dagilimlarda
# ppf bu olasiliklara karsilik gelen x degerlerini verir
x = np.linspace(st.weibull_min.ppf(0.01, k), st.weibull_min.ppf(0.99, k), 500)
pdf = st.weibull_min.pdf(x, k)
plt.plot(x, pdf, 'r-', lw=4)
# Bu dagiilma uyan bir ornek olusturalim
wb_ornek = st.weibull_min.rvs(k, size=1000)
plt.hist(wb_ornek, density=True, histtype='stepfilled', alpha=0.6)
plt.show()
# Hangi olasiliklarin hangi x degerlerine
# karsilik geldiklerine bakalim
k = 0.85
xler = st.weibull_min.ppf([0.001, 0.5, 0.999], k)
print("x: ", xler)
Başarısızlık hızının ne kadar hızlı düştüğü hem graikten hem de %1 ve %99 olasılığa karşılık gelen x değerleri arasındaki büyük farkten görülmektedir.
Poisson Dağılımı¶
Poisson rasgele değişkeni nadir gerçekleşen bir olayın bir zaman aralığında kaç kez meydana geldiğini modellemek için kullanılır. Örneğin, gökbilimdeki CCD gözlemlerinde kısa poz süreleri dahilinde dedektör üzerinde toplanan foton sayıları ya da bir web sitesinde belirli aralıklarla ziyaret edilen kullanıcı sayısı bir Poisson süreci olarak düşünülebilir.
Bir olayın bir Poisson süreci olarak tanımlanması ve bu olayın n = 1, 2, 3, ... kez gerçekleşme olasılıklarını belirlemek için Poisson dağılımının kullanılması için 3 gerekli koşul bulunmaktadır.
- Ölçümler nadir görülen olayların sayısıdır.
- Tüm ölçümler birbirinden bağımsızdır.
- Ölçümlerin görülme sıklığı ilgili zaman aralığında değişim göstermemektedir.
Poisson dağılımı olayların meydana gelme oranı ($\mu$) cinsinden tanımlanır. Bir olay bir aralıkta $n = 0, 1, 2, … $ kez oluşabilir. Bir aralıktaki ortalama olay sayısı $\lambda$ olarak belirtilir. Lambda, oran parametresi (ing. rate parameter) olarak da adlandırılan olayın gerçekleşme hızıdır. k tane olayı belirli bir aralıkta gözleme olasılığı;
$$ P(k) = e^{-\lambda}\frac{\lambda^k}{k!} $$Aslında Poisson dağılımı, Binom dağılımının, nadir gözlenen bir olayın (örneğin foton sayımı) olasılığının çok küçük olması ($p \lt\lt 1$) nedeniyle ortalama gerçekleşme sayısının ($\mu$) olası tüm olayların sayısından ($n$) çok daha küçük olduğu ($k = \mu \lt\lt n$) için özel bir durumudur. Poisson dağılımı ($p \lt\lt 1$) için Binom dağılımının bir yaklaşımı olduğundan ortalama etrafında asimetrik dağılan bir dağılımdır ve şekli 10 kez zar atma deneyinde x kez 1 gelme olasılığının dağılımına benzer.
Poisson dağılımının ortalama değeri aslında tanımında vardır. Ancak ortalama değer, $x$ için beklenen değer ($\bar{x}$) tanımından da çıkarılabilir. Standart sapması da ona eşit olduğundan Poisson dağılımı tek bir parametreyle karakterize edilir ($\mu$).
$\lambda \rightarrow \infty$ limit durumunda Poisson dağılımı Normal dağılıma yakınsar.
Bu parametre Poisson sürecine uyan bir olayın belirli bir aralıkta kaç kez gerçekleştiğni tanımlayan oran parametresidir..
scipy.stats paketi fonksiyonlarından poisson.rvs() kullanılarak Poisson dağılımına uyan rastgele bir dağılım oluşturulabilir. Olayın gerçekleştiği aralık dahilinde $\mu$ ortalama olarak olayın kaç kez gözlendiğini verir. Bu nedenle Poisson dağılımı tek parametreleri bir dağılımdır.
Örnek: Diyelim ki bir yıldız gözleminde CCD'nin herhangi bir pikseline belirli bir dalgaboyu aralığında yıldızdan saniyede ortalam 5 foton geliyor olsun. Dedektörün bu dalgaboyu aralığındaki kuantum etkinliği de %80 olsun. Yani gelen 5 fotonun 4'ü dedektör tarafından algılansın ($\mu = 4$). Gece boyunca alınan 1'er saniyelik poz süreleriyle alınan 1000 görüntüde algılanan foton sayılarının değişimi bir Poisson dağılımı sergiler.
from scipy.stats import poisson
ornek_poisson = poisson.rvs(mu=4, size=1000)
fig,ax = plt.subplots()
ax.hist(ornek_poisson)
ax.set(xlabel='Poisson Dagilimi', ylabel='Frekans')
plt.show()
Foton sayıları için gözlenme sayıları (frekans) yerine, normalize yoğunluk (o aralıktaki ölçüm sayısının toplama bölümü) verilmek istenirse pyplot.hist fonksiyonunda density parametresi $True$ değerine ayarlanır. Poisson dağılımının bir süreksiz fonksiyon olduğunu düşünerek olasılık yoğunluk fonksiyonu (PDF) yerine, istenen her değer için Poisson olasılık kütle fonksiyonunu (PMF) örnek dağılımın üzerine çizdirmek bu şekilde mümkündür. Burada PMF, fonksiyon süreksiz olduğunu vurgulamak üzere, özellikle az sayıda (10) nokta için, 'o' sembolüyle verilmiştir.
from scipy.stats import poisson
ort = 4. # oran parametresi
ornek_poisson = poisson.rvs(mu=ort, size=1000)
fig,ax = plt.subplots()
ax.hist(ornek_poisson, density = True)
# %1 ile %99 olasiliklar arasina karsilik gelen x degerleri
# arasindaki 10 nokta icin olasilik kutle fonksiyonu (pmf)
x = np.linspace(poisson.ppf(0.01, mu=ort), poisson.ppf(0.99, mu=ort),10)
pmf = poisson.pmf(x, mu=ort)
ax.plot(x, pmf, 'ro')
ax.set(xlabel='Poisson Dagilimi', ylabel='PMF')
plt.show()
Diyelim ki poz süresini uzatıp 5 saniyeye çıkarıyoruz. Ancak birim zamanda algılanan foton sayısı da değişmiyor. Buna karşın artık $\mu = 5~x~4 = 20$ foton algılarız. Doğal olarak gece boyunca alacağımız görüntü sayısı düşecektir. Poz süresini 5 kat arttırdığımız için toplam görüntü sayısı da 5 kat azalmış ve 200'e düşmüş olsun.
from scipy.stats import poisson
ornek_poisson = poisson.rvs(mu=20, size=200)
fig,ax = plt.subplots()
ax.hist(ornek_poisson)
ax.set(xlabel='Poisson Dagilimi', ylabel='Frekans')
plt.show()
Görüldüğü gibi dağılım normal dağılıma daha fazla yaklaşmıştır. Astronomide gözlemlerinde zaman çözünürlüğünden ve toplam görüntü sayısından ödün verilebilecek durumlarda poz süresinin arttırılmasıyla Poisson dağılımı gösteren foton istatistiğinin Gaussyen dağılımı yakınsaması amaçlanır. Zira görüntü işleme rutinleri, sayımların Normal dağıldığı varsayımı üzerine kuruludur.
Gamma Dağılımı¶
Gamma dağılımı, iki parametreli, sürekli olasılık dağılımlarının ailesidir. Nadiren ham formunda kullanılırken, üstel (ing. exponential), ki-kare ($\chi^2$) gibi diğer yaygın olaraak kullanılan dağılımlar ve Erlang dağılımları Gamma dağılımının özel durumlarıdır. Gamma dağılımının parametreleri, şekil parametresi, $\alpha = k$ ve oran parametresi olarak adlandırılan, ters ölçek parametresi $\beta = 1 / \theta$ 'dır. Gamma dağılımı için $\Gamma $ sembolü kullanılır ve $\Gamma (n) = (n-1)!$ olarak tanımlanır.
$$ f(x, \alpha, \beta) = \frac{\beta^{\alpha} x^{\alpha - 1} e^{-\beta x}}{\Gamma(\alpha)} $$scipy.stats paketi fonksiyonlarından gamma.rvs() kullanılarak gamma dağılımına uyan rastgele bir dağılım oluşturulabilir. $\alpha$ bir tamsayı olduğunda, gamma dağılımı Erlang dağılımına, a = 1 ise üstel dağılıma karşılık gelir.
numpy.random.gamma fonksiyonu bir gamma dağılımı oluşturmak üzere kullanılabilir. Ölçek parametresinin yukarıdaki tanımlar itibarıyla ters ölçek ($\beta$) değil, ($\beta = 1 / \theta$) tanımındaki $\theta$ parametresi olduğu dikkate alınmalıdır.
alpha = 5. # sekil parametresi
theta = 0.5 # olcek parametresi
ornek_gamma = np.random.gamma(shape=alpha, scale=theta, size=1000)
fig,ax = plt.subplots()
ax.hist(ornek_gamma, bins=50)
ax.set(xlabel='Gamma Dagilimi', ylabel='Frekans')
plt.show()
Aynı şekilde olasılık yoğunluk fonksiyonu (PDF), birikimli olasılık fonksiyonu (CDF), herhangi bir olasılığın karşılık geldiği aralık gibi gamma sürekli dağılımına ilişkin parametreler scipy.stats.gamma fonksiyonlarıyla elde edilebilir.
alpha = 3.10 # sekil parametresi
# bir gamma dagilimi olusturalim
gamma_ornek = st.gamma.rvs(alpha, size=1000)
# histogramina bakalim
plt.hist(gamma_ornek, density=True, alpha=0.8)
# bir gamma dagiliminin istatistiksel parametrelerine bakalim
ortalama, varyans, carpiklik, kurtosis = st.gamma.stats(alpha, moments='mvsk')
print("Dagilimin ortalama degeri: {:.2f}, standart sapmasi: {:.2f}, carpikligi: {:.2f} ve basikligi: {:.2f}".\
format(ortalama, varyans, carpiklik, kurtosis))
# Olasilik yogunluk fonksiyonu
# %1 ile %99 oalsilik arasini alalim. Tum dagilimlarda
# ppf bu olasiliklara karsilik gelen x degerlerini verir
x = np.linspace(st.gamma.ppf(0.01, alpha), st.gamma.ppf(0.99, alpha), 500)
pdf = st.gamma.pdf(x, alpha)
plt.xlabel("x")
plt.ylabel("PDF")
plt.plot(x, pdf, 'r-', lw=4)
# Bu dagiilma uyan bir ornek olusturalim
plt.show()
Üstel Dağılım¶
Üstel dağılım (ing. exponential distribution), bir Poisson sürecindeki (sürekli ve bağımsız ("hafızasız") bir şekilde sabit bir ortalama oranında meydana gelen) olaylar arasındaki sürenin olasılık dağılım fonksiyonu olarak tanımlanır. gamma dağılımının özel bir durumu, süreksiz bir olasılık dağılım fonksiyonu olan geometrik dağılımın sürekli durumudur. Normal, binom, Poisson, gamma gibi dağılımların üyesi olduğu üstel dağılımlar ailesinin bir üyesidir. Oran (ing. rate) parametresi olarak adlandırılan bir $\lambda$ parametresine sahiptir ve bu parametre dağılımın azalma (artma) hızını belirler. Olasılık yoğunluk fonksiyonu,
$$ p(x)~dx = \lambda~e^{-\lambda x}~dx $$ile verilir.
Burada $\lambda$ paramatresi Poisson dağılımındaki oran parametresi ($\lambda#)'ya benzer. Örneğin, CCD dedektörün bir pikseline saniyede gelen foton sayısı, Kepler Uzay Teleskobu'nun baktığı alandaki Güneş-benzeri yıldızların sayısı, bir saat içinde bir mağazaya gelen müşteri sayısı, yıllık deprem sayısı gibi nicelikler, Poisson dağılımının parametresi olan oran parametresine karşılık gelir.
Ancak, günlük hayatta olaylar arasındaki geçen süre söz konusu olduğunda, oran yerine zaman birimiyle konuşmyı yeğleriz. Örneğin, 10000 ADU sayım üretecek kadar fotonun algılanması için 20 saniye poz süresinin gerekmesi (20 / 10000 = 0.002 bir oran yerine), bir makinenin hata vermeden çalışacağı yıl sayısının 10 yıl olması (0.1 başarısızlık / yıl demek yerine) gibi. Örnekler çoğaltılabilir, bu nedenle üstel dağılımın “ortalama” terminolojisi Poisson süreçlerdeki oran parametresi $\lambda$ yerine, bunun tersi olan zaman parametresi $1 / \lambda$ üzerine kuruludur. Örneğimizden hareketle 10000 ADU 20 saniyede algılanıyorsa, "Poisson oranı" 0.002'dir, 10000 ADU algılamak için 20 saniye poz süresi vermek gerekir. Dolayısı ile $X \sim exp(0.002)$ verildiğinde kastedilen bu Poisson oranıdır.
Üstel dağılım, bir Poisson olayının gerçekleştiği iki durum arasındaki süre üzeirnden tanımlandığından bu iki olay arasındaki bekleme süresinde herhangi bir olay gerçekleşmez. Bu durum aşağıdaki şekilde ifade edilebilir.
$$ P(x = k) = \frac{\lambda^k~e^{-\lambda}}{k!} $$$$ P(x = 0) = \frac{\lambda^0~e^{-\lambda}}{0!} = e^{-\lambda} $$Tek bir zaman biriminde değil de bir $t$ zaman aralığında hiçbir olayın gerçekleşmeme olasılığı bilinmek isteniyorsa P(x=0, 1/t aralıkta) P(x=0, 1/t aralıkta) ... * P(x=0, 1/t aralıkta) çarpımı t kez yapılmalıdır. Bu durumda
$$ P(x = 0, t) = e^{-\lambda~t} $$elde edilir. Bu şekilde elde ettiğimiz olasılık "bir olay gözlendikten sonra bir sonraki olaya kadar olan bekleme süresinin t zaman biriminden büyük olma olasılığı" olup $P(t > T)$ birikimli olasılığıdır. Birikimli olasılık $P(t < T$) ile ifade olunduğundan, bu nicelik $1 - CDF$'e eşittir.
$$ P(t \le T) = 1 - P(t \gt T) = CDF = 1 - e^{-\lambda~t} $$Birikimli yoğunluk fonksiyonunun türevi olasılık yoğunluk fonksiyonunu verecektir.
$$ PDF = \frac{d~CDF}{dt} = \frac{d(1 - e^{-\lambda~t})}{dt} \rightarrow p(t)~dt = \lambda~e^{-\lambda~t}~dt $$olasılık yğunluk fonksiyonu elde edilmiş olur.
Sonuç olarak, birim zaman başına bir Poisson sürecine uyan olay sayısı bir Poisson dağılımını izliyorsa, olaylar arasındaki süre üstel dağılımı takip eder.
lamda = 1.5 # oran (olcek) parametresi
ornek_ustel = np.random.exponential(scale=lamda, size=1000)
fig,ax = plt.subplots()
ax.hist(ornek_ustel, bins=50)
ax.set(xlabel='Ustel Dagilim', ylabel='Frekans')
plt.show()
Aynı şekilde üstel bir dağılım için olasılık yoğunluk fonksiyonu (PDF), birikimli yoğunluk fonksiyonu (CDF), herhangi bir olasılığın karşılık geldiği aralık gibi sürekli dağılıma ilişkin parametreler scipy.stats.expon fonksiyonlarıyla elde edilebilir.
lamda = -1.50 # sekil parametresi
# bir ustel dagilimi olusturalim
ustel_ornek = st.expon.rvs(lamda, size=1000)
# histogramina bakalim
plt.hist(ustel_ornek, density=True, alpha=0.5)
# bir ustel dagilimin istatistiksel parametrelerine bakalim
ortalama, varyans, carpiklik, kurtosis = st.expon.stats(lamda, moments='mvsk')
print("Dagilimin ortalama degeri: {:.2f}, standart sapmasi: {:.2f}, carpikligi: {:.2f} ve basikligi: {:.2f}".\
format(ortalama, varyans, carpiklik, kurtosis))
# Olasilik yogunluk fonksiyonu
# %1 ile %99 oalsilik arasini alalim. Tum dagilimlarda
# ppf bu olasiliklara karsilik gelen x degerlerini verir
x = np.linspace(st.expon.ppf(0.01, lamda), st.expon.ppf(0.99, lamda), 500)
pdf = st.expon.pdf(x, lamda)
plt.xlabel("x")
plt.ylabel("PDF")
plt.plot(x, pdf, 'r-', lw=4)
plt.show()
Kuvvet Yasası Dağılımları¶
Kuvvet yasası (zaman zaman "ölçekleme yasası" olarak da adlandırılır), bir nicelikteki göreli bir değişikliğin, ona kuvvet yasasıyla bağlı diğer bir nicelikte orantılı bir değişme yol açtığını belirtir. Yasanın en basit örneği örneklerinden biri ters-kare yasasıdır. İki cisim arasındaki mesafe arttıkça aralarındaki çekim kuvveti mesafenin karesiyle azalır ya da iki özdeş yıldızdan yakın olanın parlaklığı, uzak olanın parlaklığından aralarındaki mesafenin karesi kadar fazladır. Yasaya uyan bir dağılımın olasılık yoğunluk fonksiyonu,
$$ p(x)~dx = a~x^{a - 1} $$ile verilir.
Örneğin $1 M_{\odot}$ ile $10 M_{\odot}$ arasındaki anakol yıldızlarının kütlelerinin olasılık dağılımı $f(M)~\alpha~M^{-3.5}$ ile bağıntılıdır, zira bir yıldızın anakol yaşı kütlesinin -3.5. kuvvetiyle orantılıdır. Yaşamakta olduğumuz CoVid19 hastalığına sebep olan SARS-CoV-2 virüsünün başlangıçtaki yayılma hızının üstel olduğu, sonrasında kuvvet yasasına uyduğu ve yine üstel olarak azaldığı tespit edilmiştir.
numpy.random.power fonksiyonu $[0, 1]$ aralığında verilen $a$ parametresinin 1 eksiğiyle (($a - 1$). kuvvetiyle) kuvvet yasasına uyan bir dağılım oluşturmak için kullanılır. $a \gt 0$ olması gerektiğine, bu anlamda fonksiyonun sadece pozitif kare yasaları (artış) için çalıştığına dikkat edilmelidir. Dağılımı ölçeklendirmek için istenen sayıyla çarpma/bölme yapılabilir.
a = 4 # dorduncu kuvvete
ornek_4kuvvet = np.random.power(a = a, size=1000)*2.5
fig,ax = plt.subplots()
ax.hist(ornek_4kuvvet, bins=50)
ax.set(xlabel='x', ylabel='Frekans')
plt.show()
Diğer dağılımlardakine benzer şekilde olasılık yoğunluk fonksiyonu (PDF), birikimli olasılık fonksiyonu (CDF), herhangi bir olasılığın karşılık geldiği aralık gibi parametreler scipy.stats.powerlaw fonksiyonlarıyla kuvvet yasasına uyan dağılımlar için de elde edilebilir.
a = 3.00 # kuvvet
# bir kuvvet yasasi dagilimi olusturalim
ornek_kuvvet = st.powerlaw.rvs(a, size=1000)
# histogramina bakalim
plt.hist(ornek_kuvvet, density=True, alpha=0.5)
# bir kuvvet yasasi dagiliminin istatistiksel parametrelerine bakalim
ortalama, varyans, carpiklik, kurtosis = st.powerlaw.stats(a = a, moments='mvsk')
print("Dagilimin ortalama degeri: {:.2f}, standart sapmasi: {:.2f}, carpikligi: {:.2f} ve basikligi: {:.2f}".\
format(ortalama, varyans, carpiklik, kurtosis))
# Olasilik yogunluk fonksiyonu
# %0.1 ile %99.9 oalsilik arasini alalim. Tum dagilimlarda
# ppf bu olasiliklara karsilik gelen x degerlerini verir
x = np.linspace(st.powerlaw.ppf(0.001, a), st.powerlaw.ppf(0.999, a), 500)
pdf = st.powerlaw.pdf(x, a)
plt.xlabel("x")
plt.ylabel("PDF")
plt.plot(x, pdf, 'r-', lw=4)
plt.show()
Eksenler logaritmik ölçekte verildiğinde logaritmanın özelliği gereği dağılımın logaritmik ölçekte eğimi kuvvete eşit şekilde doğrusal olduğu görülebilir.
plt.hist(ornek_kuvvet, density=True, alpha=0.5)
plt.plot(x, pdf, 'r-', lw=4)
plt.xscale('log')
plt.yscale('log')
plt.xlabel("x")
plt.ylabel("log(PDF)")
plt.show()
Pareto Dağılımı¶
Pareto dağılımı, kuvvet yasası dağılımları olarak bilinen ve bir $x$ değişkeninin herhangi bir değeri alma olasılığının kuvvet yasasıyla (Pareto için $a = log_45 \sim 1.16$) değiştiği dağılımlar ailesinin bir üyesidir. Adını, İtalyan mühendis ve ekonomist Vilfredo Pareto tarafından öne sürülen ve herhangi bir niceliğe etkilerin %80'inin nedenlerin sadece %20'sinden geldiğini belirten, bilimsel olmayan bir “yasa” olan Pareto İlkesi'nden (ya da Matthew İlkesi) alır. Bu ilkeye göre örneğin bir toplumunun refahının %80'i nüfusun %20'sinin elindedir.
Olasılık Yoğunluk Fonskiyonu'nun (PDF) matematiksel ifadesi,
$$ f_x (x) = \frac{\alpha x^{\alpha}_m}{x^{\alpha} + 1}, x \ge x_m $$$$ f_x (x) = 0, x \lt x_m $$ile verilir.
numpy.random.pareto fonksiyonuya bir örnek dağılım oluşturulup, Pareto fonksiyonu ifadesiyle uyumlanabileceği gibi scipy.stats.pareto fonksiyonları da dağılımın yapısına ve olasılığa ilişkin parametreleri elde etmek üzere kullanılabilir.
import numpy as np
alfa, m = 3., 2. # sekil ve mod parametreleri
s = (np.random.pareto(alfa, 1000) + 1) * m
import matplotlib.pyplot as plt
n, x, _ = plt.hist(s, 100, density=True)
pdf = alfa*m**alfa / x**(alfa+1)
plt.plot(x, max(n)*pdf/max(pdf), linewidth=2, color='r')
plt.xlabel('x')
plt.ylabel('Normalize PDF')
plt.show()
b = 5.80 # sekil parametresi
# bir Pareto dagilimi olusturalim
ornek_pareto = st.pareto.rvs(b, size=1000)
# histogramina bakalim
plt.hist(ornek_pareto, density=True, alpha=0.5)
# bir Pareto dagiliminin istatistiksel parametrelerine bakalim
ortalama, varyans, carpiklik, kurtosis = st.pareto.stats(b = b, moments='mvsk')
print("Dagilimin ortalama degeri: {:.2f}, standart sapmasi: {:.2f}, carpikligi: {:.2f} ve basikligi: {:.2f}".\
format(ortalama, varyans, carpiklik, kurtosis))
# Olasilik yogunluk fonksiyonu
# %1 ile %99 oalsilik arasini alalim. Tum dagilimlarda
# ppf bu olasiliklara karsilik gelen x degerlerini verir
x = np.linspace(st.pareto.ppf(0.01, b), st.pareto.ppf(0.99, a), 500)
pdf = st.pareto.pdf(x, a)
plt.xlabel("x")
plt.ylabel("PDF")
plt.plot(x, pdf, 'r-', lw=4)
plt.show()
# Hangi olasiligin %80'ee karsilik geldigine bakalim
# Pareto ilkesi geregi bu deger x'in
# ilk %20'lik dilimi icerisinde gerceklesir.
b = 5.80
x = st.pareto.ppf(0.20, b)
print("x: ", x)
Pareto ilkesine benzer bir ilke Zipf yasası olarak adlandırılır ve genellikle bir olayın sıralamasına göre sıklığını ifade eder. Zipf yasası, keyfi bir kitapta en sık kullanılan kelimelerin bir listesi verildiğinde, en sık kullanılan kelimenin en sık kullanılan ikinci kelimeden iki kat daha sık göründüğünü belirtir. En basit haliyle, Zipf yasası güç yasasına eşittir. Zipf yasasına uyan dağılımlar kuvvet dağılımları ailesinden Zipf Dağılımı [Zeta Dağılımı](https://en.wikipedia.org/wiki/Zeta_distribution) olarak da adlandırılır. Python'da numpy.random.zipf ve scipy.stats.zipf bu dağılımı çalışmak için gerekli fonskiyonları sunarlar.
Binom Dağılımı¶
Başarı ya da başarısızlık, kazanç ya da kayıp, "yazı" ya da "tura" gibi sadece iki sonucun mümkün olduğu ve tüm denemeler için başarı (p) ve başarısızlık olasılığının (1-p) toplamının 1 olduğu süreksizdağılıma Binom Dağılımı denir. Deneyin tekrar sayısının 1 olması durumunda elde edilen özel dağılımı Bernoulli Dağılımı adı verilir. Bununla birlikte, sonuçların eşit derecede olası olması gerekmez ve her deneme birbirinden bağımsızdır. Binom dağılımının parametreleri n ve p'dir; burada n, toplam deneme sayısıdır; p, her denemede başarı olasılığıdır. Olasılık dağılım işlevi şu şekilde verilir:
$$ f(k, n, p) = Pr (k, n, p) = {n \choose k} p^k (1 - p)^{n-k} $$Burada $n \choose k$, n'in k'lı kombinasyonlarının sayısıdır ve
$$ {n \choose k} = \frac{n!}{k!~(n-k)!} $$ile verilir.
Python'da şekil parametreleri olarak n (deneme sayısı) ve p (başarı olasılığı) argümanları olan scipy.stats modülünün binom.rvs() yöntemini kullanarak bir binom dağılımı oluşturulabilir. Dağılımın ortalama değerini değiştirmek loc parametresi kullanılır. Kodun aynı sonuçları verecek şekilde tekrarlanabilmesi için, fonksiyona random_state argümanı eklenebilir.
from scipy.stats import binom
from matplotlib import pyplot as plt
# 10 bagimsiz deney sonucunun p = 0.8 olasilik icin
# 10000 bagimsiz denemede alinma olasiliklari
ornek_binom = binom.rvs(n=10,p=0.8,size=10000)
fig,ax = plt.subplots()
ax.hist(ornek_binom, bins=20)
ax.set(xlabel='Binom Dagilimi', ylabel='Frekans')
plt.show()
Görüldüğü üzere de 10 denemenin (örneğin bir basketbol maçında denenen faul atışlarının) hepsinin birden başarıyla sonuçlanma (basketle sonuçlanması) frekansı 10000 ayrı denemede (maçta) ~1000 kez gerçekleşirken, 10 deneminin 8 isabetle sonuçlanması en sık rastlanılan durum olup ~2800 denemede (maçta) bu durumla karşılaşılmıştır. Örnekte basketbolcuların ortalama faul yüzdelerinin %80 olarak belirtilmiş olması nedeniyle bu durum doğaldır. Bir faul atışı ya başarılı ya da başarısız olacağı için bu bir Bernoulli deneyidir. Başarı olasılığı $p = 0.8 > 0.5$ olduğu için dağılım sola doğru çarpıklığa (pozitif) sahiptir (ing. left skewed). $n \rightarrow \infty$, her bir deneme için başarı olasılığı $p \rightarrow 0$ olduğu halde $np = \lambda$ sınırlı bir değerdir. $p$ ve $q = 1-p$'nin çok küçük olduğu bu durumda dağılım normal dağılıma yakınsar.
Bernoulli Dağılımı¶
Bernoulli dağılımı, binom dağılımının tek bir denemenin yapıldığı (n = 1) özel bir hali olup, olasılık kütle fonksiyonu (ing. Probability Mass Function, PMF) şu şekilde verilir:
$$ f(k, p) = p^k q^{(1-k)} = p^k (1-p)^{(1-k)} $$scipy.stats modülünün bernoulli.rvs() fonksiyonu şekil parametresi olarak p (başarı olasılığı) verilmek suretiyle bir Bernoulli Dağılımı oluşturmak için kullanılabilir.
from scipy.stats import bernoulli
# Tek bir denemenin 10000 kez yapilmasiyla
# p = 0.6 olasilikli sonucun gerceklesme sikligi dagilimi
ornek_bernoulli = bernoulli.rvs(size=10000,p=0.6)
fig,ax = plt.subplots()
ax.hist(ornek_bernoulli, bins=20)
ax.set(xlabel='Bernoulli Dagilimi', ylabel='Frekans')
plt.show()
Doğal olarak başarı = $1.0$ sonucu $p = 0.6$ olduğu için 10000 denemede 6000'e yakın sayıda gerçekleşmiştir.
Kaynaklar¶
- Python'da Dağılımlar
- Weibull Dağılımı
- Üstel Dağılım
- Üstel Dağılım - 2
- Kuvvet Yasası Dağılımları
- Binom Dağılımı
- scipy.stats Modülü Dokümantasyonu
- numpy.random Modülü Dokümantasyonu
- Ivezic, Z., Connolly, A.J., VanderPlas, J.T., Gray A., 2014, "Statistics, Data Mining, and Machine Learning in Astronomy: A Practical Guide for the Analysis of Survey Data", Princeton University Press