使用 Matlab 的中心极限定理的 PDF 和 CDF 图

Question

我正在努力绘制哪里的 PDF 和 CDF 图

Sn=X1+X2+X3+....+Xn 使用中心极限定理，其中 n = 1； 2； 3； 4； 5个； 10个； 20; 40 我将 Xi 视为 (0,3) 之间的值的均匀连续随机变量。

Here is what i have done so far - 
close all
%different sizes of input X
%N=[1 5 10 50];
N = [1 2 3 4 5 10 20 40];

%interval (1,6) for random variables
a=0;
b=3;

%to store sum of differnet sizes of input
for i=1:length(N)
    %generates uniform random numbers in the interval
    X = a + (b-a).*rand(N(i),1);
    S=zeros(1,length(X));
    S=cumsum(X);
    cd=cdf('Uniform',S,0,3);
    plot(cd);
    hold on;
end
legend('n=1','n=2','n=3','n=4','n=5','n=10','n=20','n=40');
title('CDF PLOT')
figure;

for i=1:length(N)
%generates uniform random numbers in the interval
    X = a + (b-a).*rand(N(i),1);
    S=zeros(1,length(X));
    S=cumsum(X);
    cd=pdf('Uniform',S,0,3);
    plot(cd);
    hold on;
end
legend('n=1','n=2','n=3','n=4','n=5','n=10','n=20','n=40');
title('PDF PLOT')

我的输出与我期望的相差甚远。非常感谢任何帮助。

Answer 1

这可以通过使用rand() 和cumsum() 的矢量化来完成。

例如，下面的代码生成 40 次复制，复制 10000 个 Uniform(0,3) 分布的样本，并存储在 X 中。为了满足Central Limit Theorem (CLT) 的假设，它们是independent and identically distributed (i.i.d.)。然后cumsum() 将其转换为Sn = X1 + X2 + ... 的10000 个副本，其中第一行是n = 10000 的Sn = X1 的副本，第5 行是n 的S_5 = X1 + X2 + X3 + X4 + X5 的副本。最后一行是n S_40 的副本。

% MATLAB R2019a
% Setup
N = [1:5 10 20 40];    % values of n we are interested in
LB = 0;                % lowerbound for X ~ Uniform(LB,UB)
UB = 3;                % upperbound for X ~ Uniform(LB,UB)
n = 10000;             % Number of copies (samples) for each random variable

% Generate random variates
X = LB + (UB - LB)*rand(max(N),n);     % X ~ Uniform(LB,UB)    (i.i.d.)
Sn = cumsum(X); 

您可以从图片中看到n = 2 的情况，总和确实是三角形（0,3,6）分布。对于n = 40 的情况，总和近似为正态分布（高斯），均值为 60 (40*mean(X) = 40*1.5 = 60)。这显示了probability density function (PDF) 和cumulative distribution function (CDF) 的分布收敛。

注意：CLT 通常在分布收敛到正态分布时表示，因为它已经移动，所以均值为零。通过从Sn 中减去mean(Sn) = n*mean(X) = n*0.5*(LB+UB) 来移动结果即可完成此操作。

下面的代码不是黄金标准，但它产生了图像。

figure
s(11) = subplot(6,2,1)  % n = 1
    histogram(Sn(1,:),'Normalization','pdf')
    title(s(11),'n = 1')
s(12) = subplot(6,2,2)
    cdfplot(Sn(1,:))
    title(s(12),'n = 1') 
s(21) = subplot(6,2,3)   % n = 2
    histogram(Sn(2,:),'Normalization','pdf')
    title(s(21),'n = 2')
s(22) = subplot(6,2,4)
    cdfplot(Sn(2,:))
    title(s(22),'n = 2') 
s(31) = subplot(6,2,5)  % n = 5
    histogram(Sn(5,:),'Normalization','pdf')
    title(s(31),'n = 5')
s(32) = subplot(6,2,6)
    cdfplot(Sn(5,:))
    title(s(32),'n = 5') 
s(41) = subplot(6,2,7)  % n = 10
    histogram(Sn(10,:),'Normalization','pdf')
    title(s(41),'n = 10')
s(42) = subplot(6,2,8)
    cdfplot(Sn(10,:))
    title(s(42),'n = 10') 
s(51) = subplot(6,2,9)   % n = 20
    histogram(Sn(20,:),'Normalization','pdf')
    title(s(51),'n = 20')
s(52) = subplot(6,2,10)
    cdfplot(Sn(20,:))
    title(s(52),'n = 20') 
s(61) = subplot(6,2,11)   % n = 40
    histogram(Sn(40,:),'Normalization','pdf')
    title(s(61),'n = 40')
s(62) = subplot(6,2,12)
    cdfplot(Sn(40,:))
    title(s(62),'n = 40') 
sgtitle({'PDF (left) and CDF (right) for Sn with n \in \{1, 2, 5, 10, 20, 40\}';'note different axis scales'})

for tgt = [11:10:61 12:10:62]
    xlabel(s(tgt),'Sn')
    if rem(tgt,2) == 1
        ylabel(s(tgt),'pdf')
    else                           %  rem(tgt,2) == 0
        ylabel(s(tgt),'cdf')
    end
end

用于绘图的关键函数：来自基本 MATLAB 的 histogram() 和来自统计工具箱的 cdfplot()。请注意，这可以手动完成，而无需使用几行统计工具箱来获取 cdf，然后只需调用 plot()。

---

cmets 对Sn 的方差有些担忧。

注意Sn 的方差由(n/12)*(UB-LB)^2 给出（推导如下）。蒙特卡罗模拟显示我们的Sn 样本确实具有正确的方差；事实上，随着n 变大，它会收敛于此。只需致电var(Sn(40,:))。

% with n = 10000
var(Sn(40,:))         % var(S_40) = 30   (will vary slightly depending on random seed)
(40/12)*((UB-LB)^2)   % 29.9505            

S_40 可以看出收敛性非常好：

step = 0.01;
Domain = 40:step:80;

mu = 40*(LB+UB)/2;
sigma = sqrt((40/12)*((UB-LB)^2));

figure, hold on
histogram(Sn(40,:),'Normalization','pdf')
plot(Domain,normpdf(Domain,mu,sigma),'r-','LineWidth',1.4)
ylabel('pdf')
xlabel('S_n')

Sn 的均值和方差推导：

_{对于期望（均值），第二个等式通过期望的线性成立。第三个等式成立，因为 X_i 是同分布的。}

---

它的离散版本是posted here。