第3章概率统计实例分析及MatlAb求解第3章概率统计实例分析及MatlAb求解 (1)3.1 随机变量分布与数字特征实例及MA TLAB求解 (1)3.1.1 MATLAB实现 (1)3.1.2 相关实例求解 (2)3.2 数理统计实例分析及MATLAB求解 (4)3.1.1 MATLAB实现 (4)3.1.2 相关实例求解 (4)3.3参数估计与假设检验实例分析及MATLAB求解 (5)3.1.1 MATLAB实现 (5)3.1.2 相关实例求解 (5)3.4 方差分析实例求解 (10)3.1.1 MATLAB实现 (10)3.1.2 相关实例求解 (10)3.5 判别分析应用实例及求解 (14)3.1.1 MATLAB实现 (14)3.1.2 相关实例求解 (14)3.6 聚类分析应用实例及MATLAB求解 (16)3.1.1 MATLAB实现 (16)3.1.2 相关实例求解 (16)3.1 随机变量分布与数字特征实例及MATLAB求解3.1.1 MATLAB实现用mvnpdf和mvncdf函数可以计算二维正态分布随机变量在指定位置处的概率和累积分布函数值。

利用MATLAB统计工具箱提供函数，可以比较方便地计算随机变量的分布律（概率密度函数）、分布函数及其逆累加分布函数，见附录2-1，2-2，2-3。

MATLAB中矩阵元素求期望和方差的函数分别为mean和var，若要求整个矩阵所有元素的均方差，则要使用std2函数。

随机数生成函数：rand( )和randn( )两个函数伪随机数生成函数：A=gamrnd(a,lambda,n,m) % 生成n*m的 分布的伪随机矩阵B=raylrnd(b,n,m) %生成rayleigh的伪随机数3.1.2 相关实例求解例2-1 计算服从二维正态分布的随机变量在指定范围内的累积分布函数值并绘图。

程序：%二维正态分布的随机变量在指定范围内的累积分布函数图形 mu=[0 0];sigma=[0.25 0.3;0.3 1];%协方差阵 x=-3:0.1:3;y=-3:0.2:3;[x1,y1]=meshgrid(x,y);%将平面区域网格化取值 f=mvncdf([x1(:) y1(:)],mu,sigma);%计算累积分布函数值 F=reshape(f,numel(y),numel(x));%矩阵重塑 surf(x,y,F);caxis([min(F(:))-0.5*range(F(:)),max(F(:))]);%range(x)表示最大值与最小值的差，即极差。

axis([-3 3 -3 3 0 0.5]); xlabel('x'); ylabel('y');zlabel('Probability Density');图1 二维正太分布累积分布函数值图例2-2 设X 的概率密度为⎪⎪⎪⎩⎪⎪⎪⎨⎧<<-≤≤=其他。

0;30001500,15003000;15000,1500)(22x xx x x f ，求)(X E 。

求解程序：syms xf1=x/1500^2;f2=(3000-x)/1500^2;Ex=int(x*f1,0,1500)+int(x*f2,1500,3000)运行结果：Ex =1500Ex =1/3例2-3：绘制 =0.5,1,3,5,10 时Poisson 分布的概率密度函数与概率分布函数曲线。

代码如下：x=[0:15]'; y1=[]; y2=[]; lam1=[0.5,1,3,5,10];for i=1:length(lam1)y1=[y1,poisspdf(x,lam1(i))]; y2=[y2,poisscdf(x,lam1(i))];endplot(x,y1), figure; plot(x,y2)图2 泊松分布概率密度函数图图3 泊松分布概率分布函数3.2 数理统计实例分析及MATLAB 求解3.1.1 MATLAB 实现在MATLAB 中各种随机数可以认为是独立同分布的，即简单随机样本。

常用分布的随机数产生方法，可用分布英文名称缩写加上rnd ，例如：x=betarnd(a,b,m,n) 参数为a,b 的beta 分布； x=binornd(N,p,m,n) 参数为N,p 的二项分布；3.1.2 相关实例求解例2-4：设总体密度函数cos ,,222()0,.xx f x ππ⎧-<<⎪=⎨⎪⎩其他 试从该总体中抽取容量为1000的简单随机样本。

解 利用MATLAB 编辑窗口保存以下程序，保存为ex11.mn=1000; x=zeros(1,n); k=0; while k<na=rand*pi-pi/2; b=rand/2;if b<(cos(a)/2)k=k+1;x(k)=a;endendhist(x,-pi/2:0.2:pi/2)保存完成之后，在命令窗口执行ex11，则x被赋值，且可以得到这个容量为1000的样本的直方图。

图7 直方图3.3参数估计与假设检验实例分析及MATLAB求解3.1.1 MATLAB实现3.1.2 相关实例求解例3-5:对某型号的20辆汽车记录其5L汽油的行驶里程（公里），观测数据如下：29.8 27.6 28.3 27.9 30.1 28.7 29.9 28.0 27.9 28.728.4 27.2 29.5 28.5 28.0 30.0 29.1 29.8 29.6 26.9试估计总体的均值和方差。

求解程序：%矩法估计x1=[29.8 27.6 28.3 27.9 30.1 28.7 29.9 28.0 27.9 28.7];x2=[28.4 27.2 29.5 28.5 28.0 30.0 29.1 29.8 29.6 26.9];x=[x1 x2]';muhat=mean(x)sigma2hat=moment(x,2)%样本二阶中心矩var(x,1);运行结果：muhat = 28.6950sigma2hat =0.9185例3-6：对某型号的20辆汽车记录其5L汽油的行驶里程（公里），观测数据如下：29.8 27.6 28.3 27.9 30.1 28.7 29.9 28.0 27.9 28.728.4 27.2 29.5 28.5 28.0 30.0 29.1 29.8 29.6 26.9设行驶里程服从正态分布，试用最大似然估计法估计总体的均值和方差。

求解程序：%最大似然估计x1=[29.8 27.6 28.3 27.9 30.1 28.7 29.9 28.0 27.9 28.7];x2=[28.4 27.2 29.5 28.5 28.0 30.0 29.1 29.8 29.6 26.9];x=[x1 x2]';p=mle('norm',x);muhatmle=p(1)sigma2hatmle=p(2)^2运行结果：muhatmle =28.6950sigma2hatmle =0.9185例3-7设两台车床加工同一零件，各加工8件，长度的误差为：A：-0.12 -0.80 -0.05 -0.04 -0.01 0.05 0.07 0.21 B：-1.50 -0.80 -0.40 -0.10 0.20 0.61 0.82 1.24 求方差比的置信区间。

解：用Matlab计算如下：x=[-0.12,-0.80,-0.05,-0.04,-0.01,0.05,0.07,0.21];y=[-1.50,-0.80,-0.40,-0.10,0.20,0.61, 0.82,1.24];v1=var(x); v2=var(y);c1=finv(0.025,7,7); c2=finv(0.975,7,7);a=(v1/v2)/c2; b=(v1/v2)/c1; [a,b]计算结果为： [0.0229 0.5720]结论：方差比小于1的概率至少达到了95%，说明车床A 的精度明显高。

例3-8 下面列出的是某工厂随机选取的20只零部件的装配时间（分）： 9.8 10.4 10.6 9.6 9.7 9.9 10.9 11.1 9.6 10.2 10.3 9.6 9.9 11.2 10.6 9.8 10.5 10.1 10.5 9.7设装配时间的总体服从正态分布，标准差为0.4，是否可以认为装配时间的均值在0.05的水平下不小于10。

解：H ：10<μ vs 1H ：10≥μ程序：%正态总体的方差已知时的均值检验x1=[9.8 10.4 10.6 9.6 9.7 9.9 10.9 11.1 9.6 10.2]; x2=[10.3 9.6 9.9 11.2 10.6 9.8 10.5 10.1 10.5 9.7]; x=[x1 x2]';m=10;sigma=0.4;a=0.05;[h,sig,muci]=ztest(x,m,sigma,a,1) 运行结果： h =1sig =0.01267365933873 muci =10.05287981908398 Inf因此，在0.05的水平下，可以认为装配时间的均值不小于10。

例3-9 某种电子元件的寿命X (以小时计)服从正态分布，μ和2σ均未知。

现测得16只元件的寿命如下：159 280 101 212 224 379 179 264 222 362 168 250 149 260 485 170 问是否有理由认为元件的平均寿命大于225(小时)？ 解：H ：225≤μ vs 1H ：225>μ程序：%正态总体的方差未知时的均值检验x=[159 280 101 212 224 379 179 264 222 362 168 250 149 260 485 170]; m=225;a=0.05;[h,sig,muci]=ttest(x,m,a,1) 运行结果： h=0 sig=0.2570198.2321 Inf由于sig=0.257，因此没有充分的理由认为元件的平均寿命大于225小时。

而对于H ：225≥μ vs 1H ：225<μ程序：%正态总体的方差未知时的均值检验x=[159 280 101 212 224 379 179 264 222 362 168 250 149 260 485 170]; m=225;a=0.05;[h,sig,muci]=ttest(x,m,a,-1) 运行结果： h=0 sig=0.7430 muci =-Inf 284.7679由于sig=0.743，因此更没有充分的理由认为元件的平均寿命小于225小时。

例3-10某厂铸造车间为提高铸件的耐磨性而试制了一种镍合金铸件以取代铜合金铸件，为此，从两种铸件中各独立地抽取一个容量分别为8和9的样本，测得其硬度（一种耐磨性指标）为：镍合金 76.43 76.21 73.58 69.69 65.29 70.83 82.75 72.34 铜合金 73.66 64.27 69.34 71.37 69.77 68.12 67.27 68.07 62.61根据专业经验，硬度服从正态分布，且方差保持不变，试在显著性水平05.0=α下判断镍合金的硬度是否有明显提高。