650 likes | 799 Views
数据统计分析 初级统计及回归分析 顾世梁 2008.09. 生物统计是关于试验的设计、实施,数据的收集、整理、分析和结果推论的科学。 从事试验研究,需要对处理(措施、技术)的效应给出一个明确的结论(显著与否)。 推论是先对研究对象的总体提出一种假设 (hypothesis) ,再对该假设进行测验 (test)— 以计算在假设总体中抽得实际样本 ( 统计数 ) 的概率来判断。. 1 几种常见的分布 概率计算比较复杂,生物统计中所用的概率计算主要利用变数分布进行。. 1.1 二项总体分布 ( 0 , 1 分布)
E N D
数据统计分析初级统计及回归分析顾世梁2008.09数据统计分析初级统计及回归分析顾世梁2008.09
生物统计是关于试验的设计、实施,数据的收集、整理、分析和结果推论的科学。生物统计是关于试验的设计、实施,数据的收集、整理、分析和结果推论的科学。 从事试验研究,需要对处理(措施、技术)的效应给出一个明确的结论(显著与否)。 推论是先对研究对象的总体提出一种假设(hypothesis),再对该假设进行测验(test)—以计算在假设总体中抽得实际样本(统计数)的概率来判断。
1 几种常见的分布 概率计算比较复杂,生物统计中所用的概率计算主要利用变数分布进行。 1.1 二项总体分布 (0,1 分布) 若一个总体由0,1两种元素组成,这样的总体称0,1总体。若取1的概率为p,记为P(1)=p,则P(0)=1-p=q,p+q=1.
1.2 二项分布(binomial distribution) 二项分布是指在μ=p的二项总体中,以样本容量n进行抽样,样本总和数 k (0≤k≤n))的概率分布。
1.3 普松分布(poisson distribution) 若n很大,p很小,其np=m,二项概率分布趋于普松分布。
1.4 正态分布(normal distribution) 若p接近0.5,n很大,二项概率分布趋于正态分布。
正态分布是最重要的连续性变数的分布,原因有3:正态分布是最重要的连续性变数的分布,原因有3: 1、试验研究中很多变数(性状)服从正态分布; 2、一些间断性变数在一定条件下趋于正态分布; 3、一些变数本身不服从正态,但其统计数(如平均数)在一定条件下(样本容量增大时)趋于正态分布。 这第3点是一个很重要的性质,因为我们将来对处理效应的推断,往往是以平均数(或其它统计数)进行的。在对样本容量较大的统计数进行统计推断时,可不必考虑原变数服从何种分布,统计假设测验均可在正态分布的基础上进行。
了解一个变数(或一个统计数)服从某种分布,其目标是为了计算该变数(统计数)落在某一区间的概率。P(a≤x≤b)=?了解一个变数(或一个统计数)服从某种分布,其目标是为了计算该变数(统计数)落在某一区间的概率。P(a≤x≤b)=?
1.5 学生氏 t分布( t distribution) 标准正态离差 服从正态分布。 上述u分布在实际应用中存在问题,最主要的是无法得到σ,人们自然想到用样本标准差 s代替σ计算u值,进而计算概率(假设测验)。但经抽样试验发现,这种替代是有问题的,尤其是在小样本情况下,s 的变异度较大(而σ是常量)。它直接的效果是由此算出的值比 u 的变异度大。后经WS Gosset (1908)导出了该统计数(t)的概率密度函数 f(t)。
2 统计假设测验 2.1 概念和基本步骤 我们在试验过程中获得了一个或多个样本(统计数),其目的在于推断由此代表的总体(参数)。得出处理效应存在与否的定性结论。基本过程有4步: 1)对未知总体(参数)提出假设 H0:θ=θ0, HA:θ≠θ0; H0: μ= μ0, HA: μ≠μ0; 2)设定一个否定H0假设的小概率标准(显著水平) α( α=0.05, α=0.01); 3)计算在假设条件下比实得样本(统计数)还偏的概率p。 4)根据p与α值的大小,接受或否定H0假设。
2.2 几种常用的假设测验 指的是该统计数的标准误,亦即该统计数分布的标准差。
ttest(x, m0) ttest2(x1, x1)
2.3 假设测验的本质 1)显著性 的大小是决定统计数与假设参数间、统计数间差异显著性的主要因素。试验研究中应尽量减小统计数的标准误。一是减小试验误差(s);二是增大样本容量(n)。 2)假设测验的错误 利用概率进行测验,有些情况下会犯错误。当正确的假设被否定时,就犯了弃真错误(I型错误, α错误);当错误的假设被接受时,就犯了取伪错误(II型错误, β错误)。犯两类错误的概率不同。
3 方差分析 方差分析是将多个样本作为一个整体,将总变异分解成相应变异来源的平方和和自由度,得到各变异来源方差的数量估计,用F测验鉴别样本间的差异显著性。分三个内容: 1)分解平方和自由度,计算各变异来源的方差;其中MSe(或se)比较重要,它是测验组间效应存在与否的标准; 2)F测验, F=MSt/MSe; 3)多重比较,当F测验显著,应对处理平均数的差异显著性作进一步说明。
3.1 单向分组资料的方差分析 Data structure xij为第i个处理的第j个观察值,i=1,2,…,k, j=1,2,…,n.
3.2 两向分组资料的方差分析 Data structure xij为A因素第i个水平和B因素第j个水平组合(处理)的反应量,i=1,2,…,k;j=1,2,…,n.
3.3 系统分组资料的方差分析 Data structure xijk xijk为第i组、第j亚组、第k个反应量,i=1, 2, …, l;j=1,2,…,m;k=1, 2, …, n.
较复杂的系统分组资料还可能在亚组中继续再分成小亚组(小小亚组);每一组具有不同的亚组数(mi不全相同),每一亚组具有不完全相同的观察值数目(nij不全相同)。较复杂的系统分组资料还可能在亚组中继续再分成小亚组(小小亚组);每一组具有不同的亚组数(mi不全相同),每一亚组具有不完全相同的观察值数目(nij不全相同)。 xijk为第i 组,第j亚组,第k个(处理)的反应量,i=1, 2, …, l;j=1,2,…,mi;k=1, 2, …, nij.
3.4 单因素完全随机试验资料的分析 即单向分组资料的方差分析。 3.5 单因素随机区组试验资料的分析 即两向分组资料的方差分析。 3.6 二因素随机区组试验资料的分析 A因素有a个水平,B因素有b个水平,均衡搭配时有ab个处理;r个重复(r个区组),abr个观察值。方差分析分两步:
2)构建AB两向表,按AB因素两向分解平方和、自由度。2)构建AB两向表,按AB因素两向分解平方和、自由度。
二因素、多因素完全随机试验、随机区组试验资料的方差分析均可用anovan的命令实现。二因素、多因素完全随机试验、随机区组试验资料的方差分析均可用anovan的命令实现。 格式:anovan(x, group, model)
Anovan (多因素资料的方差分析) Anovan(x, group, model) 三因素 model=[1 2 3 4 5 6 7] (三因素方差分析编码表)
3.7 一些处理效应再分解的方差分析 1)单一自由度比较; 2)其他分解的一些实例。 Lsh.m; cg.m.
如例8.1(水稻N肥试验),5个处理(ABCDE)具有SSt=301.2,dft=4,可将其进一步分解:如例8.1(水稻N肥试验),5个处理(ABCDE)具有SSt=301.2,dft=4,可将其进一步分解: ABCD vs E df1=1, SS1=198.45;AB vs CD df2=1, SS2=72.25 A vs B df3=1, SS3=12.5; C vs D df4=1, SS4=18.0
4 回归和相关分析 4.1 一元线性回归分析 对于双变数资料的回归分析,主要有三项任务: 1)建立 Y 依 X 的量化关系,即估计回归统计数和回归方程; 2)估计离回归误差,对回归方程和回归统计数进行统计假设测验; 3)回归方程的进一步利用。
模型: 据: 对Q分别对a、b求偏导并 使其为0,得正规方程组: 解得:
回归分析是用最小二乘法(least squares method)估计回归统计数B’=(a, b),使离回归平方和(Q, RSS)最小:
实例和matlab命令集 clear; clc x=[1.58, 9.98, 9.42, 1.25, .30, 2.41, 11.01, 1.85, 6.04, 5.92] y=[180, 28, 25, 117, 165, 175, 40, 160, 120, 80] x=x(:); y=y(:); n=size(y,1); SSy=var(y)*(n-1); SSx=var(x)*(n-1); xbar=mean(x); ybar=mean(y); X=[ones(n,1),x]; A=X'*X; K=X'*y; SumX=A(1,2); SumY=K(1); SumX2=A(2,2); SumXY=K(2); SP=SumXY-SumX*SumY/n C=inv(A), B=A\K, B=C*K, B=X'*X\X'*y, b=X\y Q=y'*y-B'*K, U=SSy-Q, MSQ=Q/(n-2), syx=sqrt(MSQ) F=U/MSQ; p=1-fcdf(F,1,n-2); disp(['F=',num2str(F), ' p=',num2str(p)]) sa=syx*sqrt(C(1,1)), sb=syx*sqrt(C(2,2)) ta=b(1)/sa; pa=2*tcdf(-abs(ta),n-2); disp(['ta=',num2str(ta), ' p=',num2str(pa)]) tb=b(2)/sb; pb=2*tcdf(-abs(tb),n-2); disp(['tb=',num2str(tb), ' p=',num2str(pb)]) r=corr(x,y), r2=SP^2/SSx/SSy sr=sqrt((1-r^2)/(n-2)), tr=r/sr
当其中的自变数不显著时,应将其剔除。剔除的过程应采用逐步回归的方法,即每次剔除一个偏回归平方和最小且不显著的自变数,直至所有的自变数均显著(下同)。当其中的自变数不显著时,应将其剔除。剔除的过程应采用逐步回归的方法,即每次剔除一个偏回归平方和最小且不显著的自变数,直至所有的自变数均显著(下同)。
实例和matlab命令集 clear;clc,alpha=.05; x1=[10, 9, 10, 13, 10, 10, 8, 10, 10, 10, 10, 8, 6, 8, 9]'; x2=[23, 20, 22, 21, 22, 23, 23, 24, 20, 21, 23, 21, 23, 21, 22]'; x3=[3.6,3.6,3.7,3.7,3.6,3.5,3.3,3.4,3.4,3.4,3.9,3.5,3.2,3.7,3.6]'; x4=[113, 106,111,109,110,103,100,114,104,110,104,109,114,113,105]'; y=[15.7,14.5,17.5,22.5,15.5,16.9,8.6,17,13.7,13.4,20.3,10.2,7.4,11.6,12.3]'; x=[x1,x2,x3,x4]; load regm %x=rand(100,40);y=rand(100,1); %data=xlsread('regm'); y=data(:,end);data(:,end)=[];x=data;data=[]; %data=load('regm.csv'); y=data(:,end);data(:,end)=[];x=data;data=[]; [n,m]=size(x);SSy=var(y)*(n-1); X=[ones(n,1),x]; A=X'*X;K=X'*y;C=inv(A) b=A\K,%b=C*K,b=X'*X\X'*y,b=X\y Q=y'*y-b'*K,U=SSy-Q,MSQ=Q/(n-m-1),syx=sqrt(MSQ) Fm=U/m/MSQ; p=1-fcdf(Fm,m,n-m-1);disp(['Fm=',num2str(Fm), ' p=',num2str(p)]) Up=b.*b./diag(C);Up(1)=[]; F=Up/MSQ, pr=1-fcdf(F,1,n-m-1)
for i=1:m if i<10 tr(i,:)=char(['X',num2str(i),' ']); else tr(i,:)=char(['X',num2str(i)]); end end while max(pr)>=alpha qi=find(F==min(F)); pr=1-fcdf(min(F),1,n-m-1); if pr>=alpha disp([num2str(qi),' ',num2str(min(F)),' del ',tr(qi,:)]) tr(qi,:)=[]; X(:,qi+1)=[]; m=m-1; end A=X'*X; K=X'*y; b=X\y; Q=y'*y-b'*K; MSQ=Q/(n-m-1); C=inv(A); Up=b.*b./diag(C);Up(1)=[]; F=Up/MSQ; pr=1-fcdf(F,1,n-m-1); end
disp('Last Results:') disp(' Xi bi Upi Fi p>Fi') disp(['X0 ',num2str(b(1))]) for i=1:m disp([tr(i,:),' ',num2str(b(i+1)),' ',num2str(Up(i)),' ', num2str(F(i)),' ',num2str(pr(i))]) end disp(['Error ',num2str(n-m-1),' ',num2str(Q),' ',num2str(MSQ)]) disp(['Total ',num2str(n-1),' ' num2str(SSy)]) r2=(SSy-Q)/SSy