Elementer bir anti derivativi olmayan $e^x/(x+1)$ fonsiyonunun $[0,2]$ aralığındaki integralini yaklaşık olarak hesaplayarak bu soruyu yanıtlayacağım. Bunu da, doğruluğu "büyük sayıların güçlü yasası"na dayanan bir kod yazarak yapacacağız ama önce bu yasanın ne dediğini anımsayalım.
Diyelim ki $X$ adını verdiğimiz bir deney yapıyoruz ve bu deneyin çıktıları $O\subset \mathbb{R}$ kümesinde yer alıyor. Kolaylıklık olsun diye sonuçların sürekli değil de ayrık bir küme olduğunu varsayalım. Hem bizim yöntemimizde $X$ deneyimizin yalnızca iki farklı sonucu olabilecek. Bu deneyin teorik beklentisi $E(X)$ ile gösterilir ve şu şekilde tanımlanır:
$$E(X)=\sum_{x\in O}p(x)\cdot x$$
Bu tanım içinde geçen $p(x)$ sembolü ile anlatılmak istenen, $x$ çıktısının görülme olasılığıdır. Şimdi bu deneyi tekrar tekrar yaptığımızı düşünelim, ve $N$'inci kez bu deneyi yaptığımızda elde edilen sonucu $X(N)$ ile gösterelim. Bu deneylerin ortalama sonucunu da tanımlayalım. Nasıl mı, sonuçları alıp, ortalamalarına bakalım ve bunu $X_N$ sembolü ile gösterelim. Yani,
$X_N= \frac{X(1)+X(2)+\cdots+X(N)}{N}$
olsun. Tanımıza bakarsak $X_N$ değerine ortalama çıktı desek kimse şaşırmaz herhalde. Büyük sayıların güçlü yasası, hiç deney yapılmadan tanımlanan teorik beklenti $E(X)$ ile, deney yapılarak elde edilen ortalama çıktı $X_N$ değeri arasındaki bir ilişkiyi söyler bize. O da şudur:
Büyük sayıların güçlü yasası: $N$ büyüdükçe, ortalama çıktı $X_N$ büyük olasılıkla teorik beklenti $E(X)$'e çok yakın olur.
Bunu bir örnekle açıklayalım. Diyelim ki uyanık birisiyle barbut oynuyoruz ve her zar atışından sonra $1$ ya $2$ gelirse o bize bir lira veriyor, $3,4,5,6$'dan birisi gelirse biz ona $1$ lira veriyoruz. Bu durumda $X$ deneyimizin sonucu $\{-1,+1\}$ kümesindedir ve $1/3$ olasılıkla $+1$ sonucu gerçekleşiyordur, $2/3$ olasılıkla da $-1$ sonucu. O halde
$E(X)=\frac{2}{3}(-1)+\frac{1}{3}(+1)=-\frac{1}{3}$
eşitliğinden de görülebileceği gibi bu oyundan teorik beklentimiz, $-1/3$ liradır. Bu şu demek. Eğer biz 300.000 kere zar atarsak, büyük olasılıkla ortalama olarak $1/3$ lira kaybetmiş olacağız, yani, 300.000 oyun sonunda 100.000 lira civarında para kaybetmiş olma olasılığız çok çok yüksek olacak.
Pekala, böyle bir bilgiyle nasıl integral alabiliriz? Şimdi yukarıda tanımladığım fonksiyonun grafiğine bakalım:
Sarı ile çizdiğimiz grafik, $f(x)= e^x/(x+1)$ fonksiyonunun grafiği, mavi-turkuaz karşımı düz çizgi de, benim $[0,2]$ aralığında çizdiğim uygun bir sabit $g(x)=2.5$ fonksiyonu. Bu uygun sabit $g(x)$ fonksiyonunun şöyle bir özelliği var, $f(x)$ fonksiyonumuz $g(x)$ ile sınırlı, yani sabit fonksiyonun altında kalıyor.
Şimdi bir $X$ deneyi yapacağız ve bu deneyin çıktıları ya $0$ ya da $1$ olacak. Yöntemimiz de şu olacak. $[0,2]\times [0,2.5]$ karesinden rastgele bir nokta seçeceğiz, eğer bu nokta $f(x)$ fonksiyonunun altında kalıyorsa deneyimizin sonucuna müsbet anlamında $1$ diyeceğiz, yok eğer $f(x)$'in üstünde ise, sonuç $0$ diyeceğiz.
Soru: Bu deneyin teorik beklentisi nedir?
Yanıt: Bunun için $1$ değerini bulma olasılığımızı hesaplamamız gerek. Bir dikdörtgen var, ve bu dikdörtgen içerisinden rastgele bir nokta seçiyoruz. Eğer $f$'nin altında kalan alanda ise $1$ değerini elde ediyoruz. O halde $1$ elde etme olasılığımız şöyle buluruz:
$\text{1'in gerçekleşme olasılığı}=\text{rastgele noktanın f'in altında kalan alandan olma olasılığı}=\frac{\text{f'nin altında kalan alan}}{\text{bütün alan}}$.
Şimdi teorik beklentimizi hesaplıyabiliriz:
$E(X)=\frac{\text{f'nin altında kalan alan}}{\text{bütün alan}}\times 1+\frac{\text{f'nin üstünde kalan alan}}{\text{bütün alan}}\times 0= \frac{\text{f'nin altında kalan alan}}{\text{bütün alan}}$.
Pekala, dikdörtgen içinde rastgele bir nokta seçelim ve $X(1)$'i hesaplayalım. Ya $0$ ya da $1$ olacak. Böyle böyle devam edelim ve 500000 kere rastgele nokta seçelim ve bulduğumuz değerlerin ortalamasına bakalım: Yani
$\text{500000 deney sonundaki ortalama çıktı}=\frac{\sum_{i=1}^{500000}X(i)}{500000}$
Büyük sayıların güçlü yasası der ki, çok büyük olasılıkla
$\frac{\sum_{i=1}^{500000}X(i)}{500000}\simeq E(X)=\frac{\text{f'nin altında kalan alan}}{\text{bütün alan}}$
Artık elimizde $f$'nin altında kalan alanı yaklaşık olarak veren bir formül var (büyük olasılıkla). Şimdi bu yöntemi, Python kullanarak uygulayalım.
def yaziTuraIntegrali(f, maxF, trial = 1000, a=0, b=1):
# import math as math
import random
"""f float degerler alıp float degerler alan negatif olmayan bir fonksiyon
aboveF, a,b aralığında f'nin aldığı değerlerden büyük olan bir float
f fonksiyonunun a, b aralığındaki integralinin yaklaşık değerini verir"""
bingo = 0
butunAlan = maxF * (b-a)
for i in range(trial):
# rastgele bir nokta sec
randomX = random.uniform(a,b)
randomY = random.uniform(0,maxF)
# bu nokta fonksiyonun altinda mi kaliyor kontrol et
if randomY <= f(randomX):
# altta kaliyorsa, deney sonucu 1 olacagi icin kaydet
bingo += 1
# alan ortalama sonuc formünü kullanarak alani getir
return (bingo/trial)* butunAlan
İlk önce bu fonksiyonumuzu, integralini almayı bildiğimiz fonksiyonlar için kullanıp, test edelim. Mesela $f(x)=x^2+1$ olsun. Bu fonksiyonun anti-derivativi $x^3/3+x$ ve $[0,1]$ aralığındaki integrali de $4/3=1.\overline{3}$. Ama tek bär deneme de yapmak istemiyorum, aynı deneyi tekrar tekrar yapıp bulduğum yaklaşık integrallerin de ortalamasına bakmak istiyorum, yani yazdığım programın simülasyonunu yapacağım:
def integralSimulasyon(f, simTrial, maxF, trial = 1000, a=0, b=1):
integral = 0.0
for deney in range(simTrial):
integral += yaziTuraIntegrali(f, maxF, trial, a, b)
print('Aranan integral büyük olasılıkla şu değere çok yakın: ' + str(integral/simTrial))
Şimdi bildiğimiz örnekle devam edebiliriz:
integralSimulasyon(lambda x: x**2 +1, 1, 10, trial= 10000)
Aranan integral büyük olasılıkla şu değere çok yakın: 1.3050000000000002
Hiç fena değil, ama biraz uzağız da. Deney sayısını arttıralım, bakalım büyük sayıların güçlü yasası bizi hayal kırıklığına uğratacak mı? Bu sefer kodumuzu bir kez değil, 50 kez çalıştıralım.
integralSimulasyon(lambda x: x**2 +1, 50, 10, trial= 10000)
Aranan integral büyük olasılıkla şu değere çok yakın: 1.3388000000000002
Abartalım.
integralSimulasyon(lambda x: x**2 +1, 50, 10, trial= 1000000)
Aranan integral büyük olasılıkla şu değere çok yakın: 1.3329102
Yuvalardığımız zaman $1.333$ buluyoruz ki, hiç ama hiç fena değil. Baya yakaşık bir sonuç, hiç fena değil. Kodumuza güven temin ettiğimize göre, şimdi baştaki, elementer yöntemlerle hesaplayamayacağımız fonksiyonun integralini alalım. Fonksiyonumuzun üstteni $2.5$ ile sınırlı olduğu için $fMax= 2.5$ değeri ile çağıracağız fonksiyonumuzu:
integralSimulasyon(lambda x: math.e**x/(x+1), 50, 2.5, trial= 100000, b=2)
Aranan integral büyük olasılıkla şu değere çok yakın: 2.9577509999999996
Ne güzel değil mi, taş ve sopalalarla integral alabiliyoruz!