試用Matlab軟件編制對(duì)均勻隨機(jī)數(shù)進(jìn)行性能檢驗(yàn)的各種方

1、1 試用Matlab軟件編制對(duì)均勻隨機(jī)數(shù)進(jìn)行性能檢驗(yàn)的各種方法的檢驗(yàn)，判斷程序。解：所編程序如下通過(guò)相關(guān)系數(shù)來(lái)間接檢驗(yàn)隨機(jī)數(shù)列的獨(dú)立性。（1）用來(lái)驗(yàn)數(shù)列中自第i個(gè)數(shù)開(kāi)始每m個(gè)數(shù)之間相關(guān)性的方法%設(shè)定顯著性水平為0.05，則z=1.96z=1.96;n=900;u=rand(1,n);i=1;m=2;M=floor(n-i)/m)-1;s=0;for k=0:M s=s+u(i+k*m)*u(i+(k+1)*m);endru=s/(M+1)-0.25;si=sqrt(13*M+7)/(12*(M+1);z0=abs(ru/si);if z0<=z h=0;else h=1;endh運(yùn)行程序

2、得：h=0,即接受,不拒絕獨(dú)立性假設(shè)，說(shuō)明產(chǎn)生的隨機(jī)數(shù)獨(dú)立。（2）用來(lái)驗(yàn)數(shù)列中相鄰二數(shù)相關(guān)性的方法n=900;u=rand(1,n);s=0;for i=1:n-1 s=s+u(i)*u(i+1);ends=s+u(n)*u(1);r=mean(u);p=s-n*r2;s=0;for i=1:n s=s+u(i)*u(i);endq=s-n*r2;ru=abs(p/q);if ru<2/sqrt(n) h=0;else h=1;end運(yùn)行程序得：h=0,即接受,不拒絕獨(dú)立性假設(shè)，說(shuō)明產(chǎn)生的隨機(jī)數(shù)獨(dú)立。（3）相關(guān)系數(shù)檢驗(yàn)%設(shè)定顯著性水平為0.05，則z=1.96z=1.96;n=900;u

3、=rand(1,n);r=mean(u);m=15;for j=1:m s=0; for i=1:n-j s=s+(u(i)-r)*(u(i+j)-r); end p=s/(n-1); s=0; for i=1:n s=s+(u(i)-r)2; end q=s/(n-1); ru(j)=p/q; v(j)=abs(ru(j)*sqrt(n-j); if v(j)<=z; h=0; else h=1; end h end運(yùn)行程序得：h=0,即接受,不拒絕獨(dú)立性假設(shè)，說(shuō)明產(chǎn)生的隨機(jī)數(shù)獨(dú)立。一、 均勻性檢驗(yàn)檢驗(yàn)法%均勻性檢驗(yàn)（1）kafang檢驗(yàn)法n=5000;u=rand(1,n);m=10

4、;pi=1/m;%設(shè)顯著水平，查表得kafang0.05(9)=16.919f(1:10)=0;for i=1:5000 switch floor(u(i)*10) case 0 f(1)=f(1)+1; case 1 f(2)=f(2)+1; case 2 f(3)=f(3)+1; case 3 f(4)=f(4)+1; case 4 f(5)=f(5)+1; case 5 f(6)=f(6)+1; case 6 f(7)=f(7)+1; case 7 f(8)=f(8)+1; case 8 f(9)=f(9)+1; case 9 f(10)=f(10)+1; endend s=0; for

5、 j=1:10 s=s+(f(j)-500)2/500; end s if s<16.919 h=0; else h=1; end h運(yùn)行程序得：s=9.1160: h=0。說(shuō)明產(chǎn)生的隨機(jī)樣本符合均勻性。二、 數(shù)字特征檢驗(yàn)理想均勻分布隨機(jī)數(shù)序列U(0,1)可以用三個(gè)數(shù)字特征完整地表述其統(tǒng)計(jì)特性，即均值 方差 二階原點(diǎn)矩 為了檢驗(yàn)隨機(jī)數(shù)列的獨(dú)立性與均勻性，來(lái)檢驗(yàn)樣本均值、方差、二階原點(diǎn)中心是否與理想均勻分布隨機(jī)量地相應(yīng)參數(shù)有無(wú)變化%數(shù)字特征檢驗(yàn)n=500; u=rand(1,n); r=mean(u);%求樣本均值 v=var(u);%求樣本方差 s=0; for i=1:n s=s+u(

6、i)2; end s=s/n; s %求二階矩 if r-1/2<=0.01 a1=0; else a1=1; end a1 if v-1/12<=0.01 a2=0; else a2=1; end a2 if s-1/3<=0.01 a3=0; else a3=1; end a3運(yùn)行程序得 s=0.3281;a1=0;a2=0;a3=0。說(shuō)明產(chǎn)生樣本的均值、方差、二階原點(diǎn)中心矩與理想均勻分布隨機(jī)量地相應(yīng)參數(shù)均無(wú)顯著差異2根據(jù)一單位正方形與其兩相鄰邊為半徑的1/4圓，用MCM估計(jì)值，給出求解步驟及相關(guān)程序。解：求解步驟如下： 設(shè)定所需精度； 產(chǎn)生兩個(gè) (0,1)均勻分布的隨機(jī)

7、數(shù)x和y； 判斷，若<1,則次數(shù)i+1，否則返回進(jìn)行下一次實(shí)驗(yàn)； 計(jì)算本次實(shí)驗(yàn)后的估計(jì)值pi=4*i/total,total為已進(jìn)行實(shí)驗(yàn)次數(shù)； 判斷精度，若pi和真值的差小于所設(shè)精度，停止實(shí)驗(yàn)，輸出pi,否則返回進(jìn)行下一次實(shí)驗(yàn)。3用組合法產(chǎn)生梯形分布密度的隨機(jī)變量：解：可將梯形密度函數(shù)f(x)劃分為一個(gè)梯形和一個(gè)三角形密度函數(shù)的組合。求解步驟如下： 產(chǎn)生(0,1)均勻分布的隨機(jī)數(shù)； 若<=a,產(chǎn)生一個(gè)獨(dú)立的(0,1)均勻分布的隨機(jī)數(shù)，取x= ； 若>a，產(chǎn)生一個(gè)獨(dú)立的(0,1)均勻分布的隨機(jī)數(shù)，取x=；源程序見(jiàn)附錄程序3。4試用產(chǎn)生均勻隨機(jī)數(shù)的方法產(chǎn)生2000個(gè)隨機(jī)數(shù)，通過(guò)數(shù)

8、字特征檢驗(yàn)其均勻性。解：使用MATLAB中的randn函數(shù)產(chǎn)生(0,1)均布的隨機(jī)數(shù)。檢驗(yàn)他們的均勻性只要計(jì)算出，若<,<，則認(rèn)為這些數(shù)據(jù)是均勻的。 對(duì)于本題的情況，絕大多數(shù)時(shí)候都是均勻的。源程序見(jiàn)附錄程序4。5產(chǎn)生常用連續(xù)分布隨機(jī)序列：（寫(xiě)出產(chǎn)生步驟，要求編程）1） 均勻隨機(jī)數(shù)2） 指數(shù)分布隨機(jī)數(shù)3） 瑞利分布隨機(jī)數(shù) 解：產(chǎn)生步驟如下：1） 產(chǎn)生(0,1)均布的隨機(jī)數(shù)；令；返回。2） 產(chǎn)生(0,1)均布的隨機(jī)數(shù)；令；返回。3） 產(chǎn)生(0,1)均布的隨機(jī)數(shù)；令；返回。6利用計(jì)算機(jī)產(chǎn)生一高斯分布的隨機(jī)序列，并檢查：1） 均值隨樣本大小N的變化規(guī)律；2） 方差隨樣本大小N的變化規(guī)律；3

9、） 相鄰樣本間的統(tǒng)計(jì)獨(dú)立性。解：1）由圖61可見(jiàn)，隨著N增大，序列的均值趨近與理論值0。2）由圖62可見(jiàn)，隨著N增大，序列的方差趨近與理論值1。3 由圖63可見(jiàn),隨著N增大,相鄰樣本間的獨(dú)立性減弱，趨近于理論值0。源程序見(jiàn)附錄程序6。圖 61圖 62圖 637要產(chǎn)生兩路高斯隨機(jī)數(shù)列，它們具有均值為零，方差為1，但需彼此相關(guān)，且相關(guān)系數(shù)r為定值，應(yīng)如何產(chǎn)生？給出產(chǎn)生方法，并編程實(shí)現(xiàn)。解：產(chǎn)生步驟如下： 產(chǎn)生相互獨(dú)立的(0,1)均布的隨機(jī)數(shù),； 在給定相關(guān)系數(shù)下，將相關(guān)系數(shù)陣分解為； 令(k=1,2)，則,是兩個(gè)均值為0，方差為1，相關(guān)系數(shù)為定制r的隨機(jī)序列。 源程序見(jiàn)附錄程序78試編程產(chǎn)生聯(lián)合密

10、度函數(shù)如下式所示的二維隨機(jī)數(shù)。解：的分布函數(shù)為： 則有： 令： 由于仍為0,1區(qū)間的均勻隨機(jī)數(shù)，此聯(lián)立方程之解為： 步驟： 1）產(chǎn)生0.1均勻分布隨機(jī)數(shù)列 2） 9編制應(yīng)用舉例中線譜增強(qiáng)器仿真軟件，并進(jìn)一步研究輸出信噪比與輸入信噪比的關(guān)系。仿真結(jié)果圖1 分別為原始的余弦信號(hào)、輸入信號(hào)和處理后的輸出信號(hào)。經(jīng)過(guò)自適應(yīng)線譜增強(qiáng)器，輸出信號(hào)中高斯白噪聲被較大抑制，余弦周期信號(hào)相對(duì)明顯。圖2 分別顯示了經(jīng)過(guò)1次運(yùn)算輸出的誤差和經(jīng)過(guò)500次運(yùn)算輸出的平均誤差。圖3為輸入信號(hào)和輸出信號(hào)的頻譜圖。頻譜圖說(shuō)明經(jīng)過(guò)自適應(yīng)譜線增強(qiáng)環(huán)節(jié)，白噪聲的頻譜強(qiáng)度降低，周期信號(hào)的譜線被突出了。結(jié)果表明自適應(yīng)譜線增強(qiáng)器可以和快速

11、傅立葉變換相媲美，當(dāng)未知正弦波具有一定帶寬或受調(diào)制時(shí)，則其性能優(yōu)于經(jīng)典的譜分析儀。經(jīng)過(guò)自適應(yīng)線譜增強(qiáng)器，信噪比有明顯的提高。附錄程序1一獨(dú)立性檢（1）用來(lái)驗(yàn)數(shù)列中自第i個(gè)數(shù)開(kāi)始每m個(gè)數(shù)之間相關(guān)性的方法%設(shè)定顯著性水平為0.05，則z=1.96z=1.96;n=900;u=rand(1,n);i=1;m=2;M=floor(n-i)/m)-1;s=0;for k=0:M s=s+u(i+k*m)*u(i+(k+1)*m);endru=s/(M+1)-0.25;si=sqrt(13*M+7)/(12*(M+1);z0=abs(ru/si);if z0<=z h=0;else h=1;endh

12、運(yùn)行程序得：h=0,即接受,不拒絕獨(dú)立性假設(shè)，說(shuō)明產(chǎn)生的隨機(jī)數(shù)獨(dú)立。（2）用來(lái)驗(yàn)數(shù)列中相鄰二數(shù)相關(guān)性的方法n=900;u=rand(1,n);s=0;for i=1:n-1 s=s+u(i)*u(i+1);ends=s+u(n)*u(1);r=mean(u);p=s-n*r2;s=0;for i=1:n s=s+u(i)*u(i);endq=s-n*r2;ru=abs(p/q);if ru<2/sqrt(n) h=0;else h=1;end運(yùn)行程序得：h=0,即接受,不拒絕獨(dú)立性假設(shè)，說(shuō)明產(chǎn)生的隨機(jī)數(shù)獨(dú)立。（3）相關(guān)系數(shù)檢驗(yàn)%設(shè)定顯著性水平為0.05，則z=1.96z=1.96;n=9

13、00;u=rand(1,n);r=mean(u);m=15;for j=1:m s=0; for i=1:n-j s=s+(u(i)-r)*(u(i+j)-r); end p=s/(n-1); s=0; for i=1:n s=s+(u(i)-r)2; end q=s/(n-1); ru(j)=p/q; v(j)=abs(ru(j)*sqrt(n-j); if v(j)<=z; h=0; else h=1; end h end運(yùn)行程序得：h=0,即接受,不拒絕獨(dú)立性假設(shè)，說(shuō)明產(chǎn)生的隨機(jī)數(shù)獨(dú)立。三、 均勻性檢驗(yàn)檢驗(yàn)法%均勻性檢驗(yàn)（1）kafang檢驗(yàn)法n=5000;u=rand(1,n);

14、m=10;pi=1/m;%設(shè)顯著水平，查表得kafang0.05(9)=16.919f(1:10)=0;for i=1:5000 switch floor(u(i)*10) case 0 f(1)=f(1)+1; case 1 f(2)=f(2)+1; case 2 f(3)=f(3)+1; case 3 f(4)=f(4)+1; case 4 f(5)=f(5)+1; case 5 f(6)=f(6)+1; case 6 f(7)=f(7)+1; case 7 f(8)=f(8)+1; case 8 f(9)=f(9)+1; case 9 f(10)=f(10)+1; endend s=0;

15、 for j=1:10 s=s+(f(j)-500)2/500; end s if s<16.919 h=0; else h=1; end h運(yùn)行程序得：s=9.1160: h=0。說(shuō)明產(chǎn)生的隨機(jī)樣本符合均勻性。四、 數(shù)字特征檢驗(yàn)理想均勻分布隨機(jī)數(shù)序列U(0,1)可以用三個(gè)數(shù)字特征完整地表述其統(tǒng)計(jì)特性，即均值 方差 二階原點(diǎn)矩 為了檢驗(yàn)隨機(jī)數(shù)列的獨(dú)立性與均勻性，來(lái)檢驗(yàn)樣本均值、方差、二階原點(diǎn)中心是否與理想均勻分布隨機(jī)量地相應(yīng)參數(shù)有無(wú)變化%數(shù)字特征檢驗(yàn)n=500; u=rand(1,n); r=mean(u);%求樣本均值 v=var(u);%求樣本方差 s=0; for i=1:n s=

16、s+u(i)2; end s=s/n; s %求二階矩 if r-1/2<=0.01 a1=0; else a1=1; end a1 if v-1/12<=0.01 a2=0; else a2=1; end a2 if s-1/3<=0.01 a3=0; else a3=1; end a3運(yùn)行程序得 s=0.3281;a1=0;a2=0;a3=0。說(shuō)明產(chǎn)生樣本的均值、方差、二階原點(diǎn)中心矩與理想均勻分布隨機(jī)量地相應(yīng)參數(shù)均無(wú)顯著差異程序3程序4程序5：1） u1=rand(1,2000);a=1;b=5;x1(1:2000)=0; for i=1:2000x1(i)=a+(b-a

17、)*u1(i);end2) u2=rand(1,2000);V=2;X2(1:2000)=0; For i=1:2000 X2(i)=-log(u2(i)/v;End3) u3=rand(1,2000); q=2; x3(1:2000)=0; for i=1:2000x3(i)=q*sqrt(-2*log(u3(i);end程序6clear all;close all;i = 1;k=linspace(1,25000,1000);for n = 1:length(k)x = randn(1,n);y1(i) = mean(x);y2(i) = var(x);xx1 = 0;N = length

18、(x);for j = 1:Nif j = Nx_f = x(1);else x_f = x(j+1);endx_b = x(j);xx1 = xx1+ x_b* x_f;endru(i) = (xx1- N* mean(x)2)/ (sum(x.2)- N* mean(x)2);if abs(ru)< 2/ sqrt(N);H0(i) = 1;else H0(i) = 0;endi = i+ 1;endfigure(1);plot(k, y1);grid on;xlabel('樣 本 大 小 N');ylabel('均 值');figure(2);plo

19、t(k, y2);grid on;xlabel('樣 本 大 小 N');ylabel('方 差');figure(3);plot(k, ru);grid on;xlabel('樣 本 大 小 N');ylabel('獨(dú) 立 性');程序7程序8：clcclose allclear alln=20;u1=rand(1,n);u2=rand(1,n);for i=1:n x1(i)=(u1(i)(1/3)*(u2(i)(1/2); x2(i)=1-(u1(i)(1/3);endx1程序9：clear;N=2000; %信號(hào)長(zhǎng)度F=0

20、.5; %信號(hào)數(shù)字頻率A=2; %信號(hào)幅度K=500; %迭代次數(shù)U=0.001; %迭代步長(zhǎng)delta=7;i=0:0.1:199.9;j=0:N-1;e2=zeros(1,N-10);for num=1:K e1=zeros(1,N-10); w=0 0 0 0 0 1 0 0 0 0 0; w1=fliplr(w); s=A*cos(2*pi*F*i); x=s+ sqrt(1.5)*randn(1,N); x1=zeros(1,delta),x(1:N-delta); x2=x,zeros(1,delta); %LMS for num1=1:N-10 b=x1(num1:num1+10); y(num1)=w1*b' e=x2(num1+delta)-y(num1); w1=w1+2*U*e*b; e1(num1)=e2; end e2=e2+e1;ende2=e2/K;figure(1);subplot(311);plot(j,s);axis(0,500,-2,2);title('Signal');grid on;subplot(312);plot(j,x);axis(0,800,-5,5);title('Signal+Noise');grid