實驗一拋硬幣試驗的模擬

利用python產生一系列0和1的亂數，模擬拋硬幣試驗，驗證拋一枚質地均勻的硬幣，正面向上的頻率的穩定值為0.5，
實驗步驟
（1）生成0和1的亂數序列，將其放入串列count中；也可用函式表示，
（2）統計0和1出現的次數，將其放入a中，a[0],a[1]分別表示0和1出現的次數，
（3）畫圖展示每次實驗正面向上出現的頻率

import matplotlib.pyplot as plt
import random
def doss():
    return random.randint(0,1)#回傳亂數 1為下，0為上
indices = []  #拋了多少次
freq = []  #儲存朝上的頻率
for toss_num in range(10,10001,10):
    heads = 0  #向上的次數
    for i in range(toss_num):
        if doss() == 0:
            heads+=1
    freq.append(heads/toss_num)
    indices.append(toss_num)
plt.plot(indices,freq)
plt.show()

實驗二 派的估計及其可視化

假設有一個中心在原點的邊長為2的正方形，并畫出其的內切圓，即圓的半徑為1，圓心在原點，向正方形中隨機的投點，從而計算圓周率，
使用隨機投點法計算圓周率時，只需要統計落在圓內的點和正方形內的點的個數，兩者的商為投點落在圓內的頻率，再由幾何概型，計算投點落在圓內的概率，由將該試驗大量重復的進行，根據貝努利大數定律，得到的頻率值近似等于概率，根據幾何概型有所以，從而計算出圓周率的近似值，

• 實驗步驟：
• 1）生成兩個系列的[-1，1]的亂數，
• 2）如果點在圓內，把點標為藍色；否則標綠色，
• 3）統計藍色點的個數，計算頻率，從而算出圓周率的近似值，

import numpy as np
import matplotlib.pyplot as plt

# 畫單位圓
theta = np.linspace(0, 2 * np.pi, 200)
x = np.cos(theta)
y = np.sin(theta)
fig, ax = plt.subplots(figsize=(4, 4))
ax.plot(x, y, color="darkred")
plt.xticks([])  # 去掉坐標軸的刻度等
plt.yticks([])  # 去掉坐標軸的刻度等

N = 1000
count = 0
for i in range(1, N + 1):
    x = np.random.uniform(-1, 1)  # 產生服從[-1,1]上的均勻分布的N個亂數
    y = np.random.uniform(-1, 1)  # 產生服從[-1,1]上的均勻分布的N個亂數
    if x ** 2 + y ** 2 < 1:
        count += 1
        plt.plot(x, y, 'bo')  # 如果點在圓內，畫藍色
    else:
        plt.plot(x, y, 'go')  # 不在圓內，畫綠色

p = count / N
plt.show()
pai = 4 * p
print(pai)

實驗三 二項分布的應用

有一繁忙的汽車站，每天有大量汽車通過，設每輛車在一天的某時段出事故的概率0.0001，在某天的該時段內有1000輛汽車通過，問出事故的次數不小于2的概率是多少（利用泊松定理）？

• 1）先使用二項分布的分布函式計算X>1的概率，
• 2）再使用公式計算X>1的概率
• 對比兩個概率值，
• 3）使用泊松分布的分布函式計算X>1的概率，
• 4）使用公式計算X>1的概率，
• 5）對比使用兩個分布計算出來概率值，

import numpy as np
from scipy.stats import binom           #匯入binom
from matplotlib import pyplot as plt

# numpy.random.binomial(n,p,size=None)

# 3個引數：n表示伯努利試驗次數，p表示伯努利試驗得到正例的概率，size表示采樣次數；回傳結果為出現正例的次數k，

#cdf(k, n, p) k：分布函式自變數，n，p：分布引數  累積概率函式（分布函式）F ( k ) 
#sf(k, n, p) k：分布函式自變數，n，p：分布引數  殘存函式 1-F（k）

prob=binom.sf(k=1, n=1000, p=0.0001)    #計算P(X<=1)
print(prob)

from scipy.stats import poisson #匯入poisson
prob1=poisson.sf(k=1, mu=1)   #計算引數為1時概率P(X>1)
print(prob1)

實驗四 正態分布的計算

設某品牌瓶裝礦泉水的標準容量是500毫升，設每瓶容量X(以毫升)是隨機變數，X~N(500，25)，求：

• (1) 隨機抽查一瓶，其容量大于510毫升的概率；
• (2) 隨機抽查一瓶，其容量與標準容量之差的絕對值在8毫升之內的概率；
• (3) 求常數C，使每瓶的容量小于C的概率為0.05，
• 實驗步驟：
• （1）用norm.cdf()計算正態分布的分布函式值，normcdf(x,mu,sigma)
  其中norm.cdf(x,mu,sigma)分別表示要求分布函式的x的值，mu為期望，sigma為標準差，
• （2）使用norm.cdf(510，500，5)計算，
• （3），使用norm.cdf(508，500，5)- norm.cdf(5492，500，5)計算，
• （4）norm.ppf(alpha,mu,sigma)計算α分位數，

from scipy.stats import norm
prob1=norm.cdf(510,500,5)
print(prob1)

prob2=norm.cdf(508,500,5)- norm.cdf(492,500,5)
print(prob2)

alpha=norm.ppf(0.05,500,5)
print(alpha)

實驗五 二維連續型分布的相關系數

設隨機變數聯合密度函式為[外鏈圖片轉存失敗,源站可能有防盜鏈機制,建議將圖片保存下來直接上傳(img-ZvXAyQNs-1639205559720)(attachment:image.gif)]，
求X、Y的期望、相關系數，
實驗步驟：

• （1）設5個函式y1,y2,y3,y4,y5，分別表示!
• （2）令y1中的兩個變數分別從0到2積分，使用from scipy.integrate import dblquad匯入二重積分的包，而計算二重積分的格式為dblquad(y1,0,2,lambda g:0,lambda h:2)，其中前面的0，2表示x的積分限，lambda g:0,lambda h:2表示y的積分限，而且dblquad回傳兩個值，一個是積分值，另一個是誤差，可設兩個變數接收結果，后面需要使用時僅用積分值即可，最終得到E(X)，
• （3）令y2中的兩個變數分別從0到2積分，得到E(Y)，令y3中的兩個變數分別從0到2積分，得到E(X2)，令y4中的兩個變數分別從0到2積分，得到E(Y2)，令y5中的兩個變數分別從0到2積分，得到E(XY)，
• （4）使用方差的計算公式，計算D(X)，D(Y)，使用協方差的計算公式，計算協方差，再使用，求出相關系數，

from scipy.integrate import dblquad
import math

Ex, err1 = dblquad(lambda y,x:1/8*x*(x+y),0,2,lambda g:0,lambda h:2)
print(Ex)

Ey, err2 = dblquad(lambda y,x:1/8*y*(x+y),0,2,lambda g:0,lambda h:2)
print(Ey)

Ex2, err3 = dblquad(lambda y,x:1/8*x*x*(x+y),0,2,lambda g:0,lambda h:2)
print(Ex2)

Ey2, err4 = dblquad(lambda y,x:1/8*y*y*(x+y),0,2,lambda g:0,lambda h:2)
print(Ey2)

Exy, err3 = dblquad(lambda y,x:1/8*x*y*(x+y),0,2,lambda g:0,lambda h:2)
print(Exy)

Dx = (Ex2-Ex)**2
print(Dx)

Dy = (Ey2-Ey)**2
print(Dy)

cov = Exy-Ex*Ey
print(cov)

p = cov/math.sqrt(Dx*Dy)
print(p)

本文用到的文章
統計函式庫scipy
二重積分求法1
二重積分求法2
常用python模擬分布及其可視化1
常用python模擬分布及其可視化2
🌟離散型隨機變數（二項分布和泊松分布及其應用）
蒙特卡洛方法計算圓周率并將撒點分布和π估計值可視化