Mathematica 教學：Central Limit Theorem

Clear["Global`*"];
SetOptions[Plot,{Frame->True,Axes->None,PlotStyle->Red}];
SetOptions[Histogram,{Frame->True,Axes->None,ImageSize->350}];

mySumStdDist[dist_,parameters__,size_,repeat_,scale_]:=
Block[{x,y,domain,mydist,mysumdist,mypdf,normaldist,mu,sigma,
range,data1,data2,data3,data4,data1std,hist1,hist2,hist3,
hist4,normalplot,normalstdplot,histscale},
mydist=Apply[dist,parameters];

mysumdist=TransformedDistribution[Total[y[#]&/@Range[size]],
(y[#]\[Distributed]mydist)&/@Range[size]];

mu=Mean[mysumdist];
sigma=StandardDeviation[mysumdist];
domain=DistributionDomain[mysumdist];
range=If[Head[domain]===Interval,domain[[1]],domain];
range={Max[range[[1]],mu-3sigma],Min[range[[-1]],mu+3sigma]};

data1=RandomVariate[mysumdist,size];
data1std=Standardize@data1;
data2=RandomVariate[mydist,{repeat,size}];
data3=Total/@data2;
data4=Standardize[Mean/@data2];

(*直方圖*)

hist1=Histogram[data1,{Min[data1],Max[data1],
(Max[data1]-Min[data1])/20.},"PDF"];

hist2=Histogram[data1std,{Min[data1std],Max[data1std],
(Max[data1std]-Min[data1std])/20.},"PDF"];

hist3=Histogram[data3,{Min[data3],Max[data3],
(Max[data3]-Min[data3])/20.},"PDF"];

hist4=Histogram[data4,{Min[data4],Max[data4],
(Max[data4]-Min[data4])/20.},"PDF"];

normalplot=Plot[PDF[NormalDistribution[mu,sigma],x],{x,Min[data1],Max[data1]}];

normalstdplot=Plot[PDF[NormalDistribution[0,1],x],{x,-4,4}];

(*避免直方圖機率小於設定scale*)

histscale=Max[scale,Max[Cases[#,RectangleBox[a__,b__,___]:>b[[-1]],
Infinity]]]&/@{hist1,hist2,hist3,hist4};

(*圖形輸出*)

GraphicsGrid[
{{
Show[hist1,normalplot,
PlotRange->{0,histscale[[1]]},
PlotLabel->PlotLabel->"The distribution of Sum Xi"],

Show[hist2,normalstdplot,
PlotRange->{0,histscale[[2]]},
PlotLabel->"The distribution of Standardize Sum Xi"]},

{Show[hist3,normalplot,
PlotRange->{0,histscale[[3]]},
PlotLabel->"The central limit theorem for Sum Xi"],

Show[hist4,normalstdplot,
PlotRange->{0,histscale[[4]]},
PlotLabel->"The central limit theorem for Standardize sample mean"]
}},ImageSize->700]
];


Manipulate[
mySumStdDist[BinomialDistribution,{n,p},size,repeat,scale],
{{n,10},10,100,1,Appearance->"Labeled"},
{{p,0.2},0.1,0.9,0.1,Appearance->"Labeled"},
{{size,100},30,1000,1,Appearance->"Labeled"},
{{repeat,100},10,1000,1,Appearance->"Labeled"},
{{scale,0.5},0.01,1,0.01,Appearance->"Labeled"}
]

Manipulate[
mySumStdDist[GammaDistribution,{a,b},size,repeat,scale],
{{a,10},1,100,1,Appearance->"Labeled"},
{{b,0.5},0.1,1.5,0.1,Appearance->"Labeled"},
{{size,100},30,1000,1,Appearance->"Labeled"},
{{repeat,100},10,1000,1,Appearance->"Labeled"},
{{scale,0.5},0.01,1,0.01,Appearance->"Labeled"}
]

Tags: Mathematica 教學, Mathematica 教學 畫圖, Mathematica 教學 程式, Mathematica教學講義

• 編輯文章
• 此文章網址
• 轉寄此文章
• 目前有0 意見