第10章回归分析

第10章回归分析本章主要内容概述Corr过程线性相关分析Reg过程直线回归多元线性回归岭回归Nlin过程非线性回归分析对数曲线拟合Logistic过程主成分回归概述变量之间的关系有两类：一类是变量间存在着完全确定性的关系，可以用精确的数学表达式来表示。如长方形的面积（S）与长（a）和宽（b）的关系可以表达为：S=ab。它们之间的关系是确定性的，只要知道了其中两个变量的值就可以精确地计算出另一个变量的值，这类变量间的关系称为函数关系，模型为机理模型。概述变量之间的另一类关系是不完全确定的关系，无法用精确的数学公式来表示。如儿童生长时的身高和体重、作物的产量和肥料用量、植物病害的发生程度与气候条件之间的关系。这类关系在生物界中大量存在，只能用统计学的方法研究变量间的相关关系，建立数值模拟模型。概述统计学上采用回归分析（regressionanalysis）研究相关变量间的因果关系。表示原因的变量称为自变量，表示结果的变量称为依变量。研究“一因一果”，即一个自变量与一个依变量的回归分析称为一元回归分析；研究“多因一果”，即多个自变量与一个依变量的回归分析称为多元回归分析。概述统计学上采用相关分析（correlationanalysis）研究变量之间的相关关系。对两个变量间的直线关系进行相关分析称为简单相关分析（也叫线性相关分析）；对多个变量进行相关分析时：研究一个变量与多个变量间的线性相关称为复相关分析；研究其余变量保持不变的情况下两个变量间的线性相关称为偏相关分析；研究两组变量间线性相关称为典型相关分析。概述SAS系统中，适用回归分析的过程：过程功能REG执行普通线性回归分析。适用于各式的输入输出格式并有诊断性以及简化模型的功能。PROBIT执行概率回归分析或逻辑斯谛的回归分析。这个程序所处理的数据通常含二分(或二分以上)的依变量以及数个连续的自变量。NLIN建立非线性的回归模型。概述过程功能LOGISTIC执行逻辑斯谛的回归分析，分析方式含逐步回归分析以及各式的诊断统计值。ORTHOREG使用Gentleman-Givens的计算程序来估计回归模型中的参数值。适用于估计值之标准误差较大的数据。RSREG建立二项式反应面(Response-Surface)的回归模型。GLM最普通的线性分析，自变量可以是类别变量或多项式。概述几种回归过程相同的输出数据用最小二乘法所估计的参数值(如：b0，b1....)；误差方差的估计值；参数估计值的标准误或方差；有关参数的假设(如H0：β0=0)检验；各种预测值及其误差；回归方程显著性的检验。概述相异的诊断功能REG程序提供共线性(Collinearity)的诊断，这个诊断探讨自变量间相关的程度及可能造成的影响。REG，LOGISTIC及RSREG三个程序提供影响度诊断以决定各观察体对参数估计值、误差的平方和(SSE)及预测值等的影响。LOGISTIC程序也有这种功能，不过其分析原理是采最大似然率法。PROBIT与RSREG两程序提供回归模型精确度(Accuracy)的诊断，所用的方法是比较误差的方差及其估计值。REG程序提供时间序列分析(TimeSeriesAnalysis)的诊断，特别是有关时间的误差以及误差间彼此的相关。Corr过程Corr过程存在于SAS的base模块，可以计算：Pearson积矩相关系数；Spearman秩相关系数；Kendall‘stau-b统计量；Hoeffding‘s独立性分析统计量D；Pearson、Spearman以及Kendall偏相关系数；另外，它还对用于估计可靠性的Cronbach系数α进行计算。Corr过程Corr过程格式：PROCCORR<options>;

BY<DESCENDING>var-1<...<DESCENDING>var-n><NOTSORTED>;

FREQ

variable;

PARTIAL

variables;

VAR

variables;

WEIGHT

variable;

WITHvariables;Run;Corr过程选项功能alpha计算并输出Cronbach系数αbest=n对于每个变量输出所指定的n个相关系数cov指定计算并输出协方差矩阵csscp指定计算并输出校正离均差平方和及离均差积和矩阵data=指定要进行分析的数据集exclnpwgt将weight语句指定变量值非正的记录排除出分析过程hoeffding指定计算并输出Hoeffding'sD统计量Corr过程选项功能kendall指定计算并输出Kendalltau-b系数nocorr禁止Pearson相关系数的计算和输出nomiss将含有缺失值的记录排除出计算过程之外noprob禁止相关系数对应的P值的计算和输出nosimple禁止描述性统计结果的输出outh=创建含有Hoeffding’sD统计量的输出数据集outk=创建含有Kendall相关统计量的输出数据集Outp=创建含有Pearson相关统计量的输出数据集

outs=创建含有Spearman相关统计量的输出数据集Corr过程选项功能pearson指定计算并输出Pearson积矩相关系数rank将相关系数按照其绝对值从大到小的顺序显示在输出结果中singular=p指定判断某变量为奇异的标准spearman计算并输出Spearman秩相关系数sscp计算并输出离均差平方和及离均差积和矩阵vardef=分母指定计算方差、标准差及协方差时所用的分母项，默认值为自由度dfouts=创建含有Spearman相关统计量的输出数据集线性相关分析线性相关分析：是在（X，Y）服从正态分布（双变量正态分布）的假定下，用线性相关的统计分析方法研究两变量间数量变化关系的密切程度和方向，称为线性相关分析。对样本的要求:线性：linear独立：independent正态：normal等方差：equalvariance线性相关分析线性相关分析线性相关分析线性相关的性质可由散点图直观的说明，按图中散点的分布，可归纳为以下几种情况：A0<ρ<1B-1<ρ<0Cρ=1Dρ=-1Eρ=0Fρ=0Gρ=0Hρ=0线性相关分析公式的解释：r的符号（方向）：

当X，Y同向变化时，一、三象限的点多于二、四象限的点，r为正，两变量正相关，0<r<1；当X，Y反向变化时，二、四象限的点多于一、三象限的点，r为负，两变量负相关，-1<r<0；当点在一、三象限和二、四象限均匀分布时，正负抵消，r=0，呈零相关，X、Y没有线性关系。相关系数r的大小：r的绝对值越大，说明线性相关程度越好，散点越接近一条直线；绝对值越小，说明线性相关程度越差。线性相关分析示例10001（李春喜，王文林等编著《生物统计学》114页）：有人研究了粘虫孵化历期平均温度（x，℃）与历期天数（y，天）之间的关系。试进行皮尔逊积差相关分析。X11.814.715.616.817.118.819.520.4y30.117.316.713.611.910.78.36.7线性相关分析对资料进行相关分析之前，我们可以先对其绘制散点图，以考察两变量的真实变化关系。我们可以应用plot过程或gplot过程来执行绘制散点图的功能。散点图完成后再计算变量之间的相关系数，对相关系数进行假设检验，以量化形式表示变量间的相关关系。线性相关分析/*SasProg10001.Sas*/dataD10001;

inputxy@@;cards;11.830.114.717.315.616.716.813.617.111.918.810.719.58.320.46.7;proc

gplot;

ploty*x='*';run;proc

corr;

varxy;run;Pearson相关分析是Corr过程的默认方式；用with语句可以对特定的变量对进行相关分析，可将corr过程中间的语句改为：

varx；withy；线性相关分析结果解释：首先给出两变量的描述性统计量；然后给出变量的相关系数矩阵（var语句所列变量中任两者之间的相关系数）；由于抽样误差的存在，从ρ=0的双变量总体中抽出的样本相关系数r不一定等于0。为了判断r所代表的总体是否存在直线相关，必须测定r来自ρ=0总体的概率，因此，需进行相关系数的假设检验。假设为H0：ρ=0，检验结果（仅给出P值）列在相应相关系数的下面。本例p<0.001，ρ=0是不可能的，x和y之间存在着极显著的直线相关关系。线性相关分析示例10002（西北农学院、华南农业大学主编《农业化学研究法》177页）：大豆脂肪含量（x1，%）和蛋白质含量（x2，%）的测定结果。试进行皮尔逊积差相关分析。X116.517.518.519.520.521.522.523.524.5X243.542.642.640.640.338.737.236.034.0线性相关分析/*SasProg10002.Sas*/dataD10002;

inputx1x2@@;cards;16.543.517.542.618.542.619.540.620.540.321.538.722.537.223.536.024.534.0;proc

gplot;

plotx2*x1='*';run;proc

corr;

varx1x2;run;线性相关分析示例10003：用60Co对狗造成急性放射病，对照射后5天时的健康状况进行综合评分，并记录其存活天数。试作Spearman秩（等级）相关分析。综合评分79809190708792存活天数45301624282514线性相关分析当两变量不符合双变量正态分布的假设时，需用Spearman秩相关来描述变量间的相互变化关系。此时，散点图上散点的分布形态不能完全描述两变量间的相关关系，故此时一般不需再绘制散点图。线性相关分析/*SasProg10003.Sas*/dataD10003;

inputxy@@;cards;7945803091169024702887259214;proc

corr

spearman

nosimple;

varxy;run;Proccorr语句的spearman选项即要求用Spearman秩相关来进行数据分析，nosimple选项则用来禁止对变量描述性统计量的输出。REG过程REG程序将参数估计值带入线性回归模型中，以便执行回归分析的预测。这些参数估计值是由最小误差平方法所导出的。REG程序是SAS所有回归分析程序中用途最广泛的一种(好比GLM程序是SAS所有方差分析程序中用途最广泛的一种)。其它回归分析的程序各有其特殊的用途。REG过程PROCREG<options>;<label:>MODELdependents=<regressors></options>;BYvariables;FREQvariable;IDvariables;VARvariables;WEIGHTvariable;ADDvariables;DELETEvariables;<label:>MTEST<equation<,...,equation>></options>;OUTPUT<OUT=dataset>keyword=names<...keyword=names>;PAINT<condition

|ALLOBS></options>|<STATUS|UNDO>;PLOT<yvar*xvar><=symbol><...yvar*xvar><=symbol></option>;PRINT<options><ANOVA><MODELDATA>;REFIT;RESTRICTequation,...,equation;REWEIGHT<condition|ALLOBS></options>|<STATUS|UNDO>;<label:>TESTequation,<,...,equation></options>;Run;REG过程Reg过程的Model语句Model语句用以指定所要拟合的回归模型。格式：<label:>MODEL

dependents=<regressors></options>;最前面的标签为可选项，可以是不超过8个字符的字符串，用来对定义的模型进行标识，以便于在结果中分辨不同的模型，一般情况下系统会以默认的方式对模型进行标识。REG过程Reg过程的Model语句关键字Model后所列的是模型表达式，和方差分析中anova过程的model语句相似。模型表达式中等号的左边为依变量，等号的右边为自变量列表，自变量间以空格相分隔。Model语句中不允许出现组合型变量。如果要用到几个变量产生的综合变量，必须在数据步完成新变量的创建过程。REG过程常用的Model语句选项：选项功能adjrsq对每一个所选择的模型计算其自由度校正的R2b对于以R2、校正R2、Cp值为参数的选择模型，计算其回归系数clb计算并输出参数估计值的可信区间上、下限cli计算并输出单个预测值的可信区间上、下限clm计算并输出每条观测下应变量期望值（均数）的可信区间上、下限collin对自变量之间的共线性进行分析REG过程选项功能collinoint对自变量之间的共线性进行分析，不包括截距项cp对每一个模型计算Mallow’sCp统计量include=n指定每次模型拟合必须包含自变量中的前n个influence针对每一条观测，分析其对参数估计和预测值的影响maxstep=n限定逐步回归分析时最多进行的步数为nnoint要求模型拟合时不包含截距项p计算每一条观测应变量的预测值并输出到结果REG过程选项功能partial对每一个自变量（包括截距项）绘制对于应变量的偏回归残差图r进行残差分析并显示在结果中sbc计算每个模型的SBC统计量并显示在结果中selection=指定模型选择的方法，可以是前进法（forward）、后退法（backward）、逐步法（stepwise）等等slentry=指定前进法和逐步法时变量进入模型的显著性水平，默认值前进法为0.5，逐步法为0.15slstay指定后退法和逐步法时变量留在模型内的显著性水平，默认值后退法为0.1，逐步法为0.15REG过程选项功能sp对每个模型计算Sp统计量ss1将I型平方和（SS1）与各参数估计值一并列出ss2将II型平方和（SS2）与各参数估计值一并列出sse计算并显示每个模型的误差平方和start=s在模型选择模式下，指定开始时包含在模型中的自变量个数（s个）stop=s指定包含在模型中的自变量个数为s个时停止模型选择过程xpx计算并显示自变量交叉积和矩阵（X'X）直线回归示例10R01（西北农学院、华南农业大学主编.农业化学研究法.北京：农业出版社.1980年12月第1版）：从土壤学理论知道，土壤中的氮素主要以有机氮的形态存在，因此，土壤全氮量与土壤有机质之间，应该有一定的联系。50种水稻土实测得到的有机质含量（%）与全氮含量（%）的记录（SasData10R01.txt），试对这批数据进行回归分析。直线回归/*SasProg10R01.Sas*/datad10R01;

infile

'e:\sas\sasdata\sasdata10R01.txt';

inputxy;proc

gplot;

ploty*x;proc

reg;

modely=x;

ploty*x;run;直线回归结果解释：Model:MODEL1：模型标签为“MODEL1”（程序未进行指定，此处为系统默认分配的标签）。DependentVariable:y：依变量为y。AnalysisofVariance：方差分析表。RootMSE：误差均方平方根。R-Square：决定系数。DependentMean：依变量均数。AdjR-Sq：校正的决定系数。直线回归结果解释：CoeffVar：变异系数。ParameterEstimates：模型的参数估计。分别给出截距项和自变量回归系数等的估计值以及对应参数为0的假设测验结果。本例回归方程为：y=0.00093116+0.05431x直线回归直线回归示例10R02（蒋明.SAS软件在回归分析中的应用.计算机与农业.2003,(8):15-16）：国营大岭农场橡胶树大型系比二组无性系，1960年刺检干胶量(mg)与1965年正式割胶产量(g)数据（SasData10004.txt)，试求正式割胶量(x)回归于刺检干胶量(y)的回归方程。直线回归/*SasProg10R02.Sas*/datad10R02;

infile

'e:\sas\sasdata\sasdata10R02.txt';

inputxy;proc

reg;

modely=x;run;本例回归方程为：y=1.96825+0.07761x多元线性回归示例10M01（蒋明.SAS软件在回归分析中的应用.计算机与农业.2003,(8):15-16）：研究橡胶树白粉病的始病期(x1)、嫩叶历期(x2)、发病期雨天数(x3)与最终病情指数(y)的关系，18年的白粉病预报要素和最终病情指数级值数据(SasData10M01.txt)，试求回归方程。多元线性回归/*SasProg10m01.Sas*/datad10m01;

infile

'e:\sas\sasdata\sasdata10m01.txt';

inputx1x2x3y;proc

reg;

modely=x1x2x3;run;多元线性回归模型选择：上例整个模型具有统计学意义（方差分析表中F检验p<0.0001），但部分自变量回归系数的检验结果显示没有统计学意义（p>0.05）。我们看到，并不是所有自变量对依变量的回归作用都具有统计学意义，这时我们需要建立一个最优模型。SAS有关模型的选择方法多达9种。多元线性回归逐步回归：/*SasProg10m02.Sas*/datad10m01;

infile

'e:\sas\sasdata\sasdata10m01.txt';

inputx1x2x3y;proc

reg;

modely=x1x2x3/selection=stepwiseslentry=0.05

slstay=0.05;run;逐步法选择模型，自变量进入模型和剔除出模型的标准都定为α=0.05。多元线性回归示例10M03：某学校20名一年级女大学生体重（x1,公斤）、胸围（x2,厘米）、肩宽（x3,厘米）及肺活量（y,升）实测值，试对影响女大学生肺活量的有关因素作逐步回归分析。数据文件SasData10M03.txt多元线性回归/*SasProg10m03.Sas*/datad10m03;

infile

'e:\sas\sasdata\sasdata10m03.txt';

inputnox1x2x3y;proc

reg;

modely=x1x2x3/selection=stepwiseslentry=0.05

slstay=0.05;run;quit;多元线性回归实例10M04（陈立，徐汉虹，赵善欢.获取复配农药最佳增效配方的一种简易方法.生物数学学报.2001,16(4):456-461）：研究氯氰菊酯和喹硫磷混配对斜纹夜蛾(Spodopteralitura)的增效作用，以共毒系数的倒数为目标函数寻找最优配方。药剂为91.6%氯氰菊酯原药(广东省中山凯达精细化工股份有限公司)和40%喹硫磷乳油(福建三明农药厂)。采用二次回归通用旋转组合设计。(1)首先测定氯氰菊酯和喹硫磷各自的毒力回归线，以Lc25的近似值为零水平。然后用以2为底的对数尺度2r编码(见下表)，杀虫剂使用浓度(Ni)与编码值(xi)之间的转换关系为多元线性回归实例10M04（陈立，徐汉虹，赵善欢.获取复配农药最佳增效配方的一种简易方法.生物数学学报.2001,16(4):456-461）：因子及编码表供试昆虫为室内饲养的斜纹夜娥2龄中后期幼虫。采用Potter塔喷雾法处理试虫，24小时后检查试虫死亡情况。编码值氯氰菊酯(x1)喹硫磷(x2)浓度(μg/ml)浓度对数浓度(μg/ml)浓度对数1.41439.971.601766.621.82361301.4771501.69900151.1761251.3979-17.50.875112.51.0969-1.4145.630.75049.380.9723多元线性回归/*SasProg10M04.Sas*/dataxz;

inputx1-x5y;cards;-1-11114.9164-11-1115.34061-1-1115.1679111116.2816-1.41400204.83251.41400205.52440-1.4140024.916401.4140025.8416000005.2533000005.4308000005.4308000005.3406000005.5244;proc

reg;

modely=x1-x5;run;计算斜纹夜蛾死亡率(y)与氯氰菊酯(x1)和喹硫磷(x2)浓度的二次回归方程。y=5.39598+0.27140x1

+0.35582x2

+0.17237x1x2

-0.07179x12

+0.02849x22多元线性回归示例10M05：预测人体吸入氧气的效率本例数据（SasData10M05.txt）来自一群中年男士的健康状态资料，由Linnerud提供。每一名男士提供七个数据，分别是：年龄（AGE），体重（WEIGHT），吸氧的效率（OXY），跑1.5英里所需的时间（以分钟计，RUNTIME），休息时的心跳（RSTPULSE），跑步时的心跳率（RUNPULSE）和最高心跳率（MAXPULSE）。其中，吸氧效率（OXY）是依变量，另外六个均是自变量。分析的过程是，先用逐步排除法，再用最大相关法，以便找出一个又精简又有效的递归模型。多元线性回归/*Sasprog10M05.Sas*/DATAD10M05;

infile

'e:\Sas\SasData\SasData10M05.txt';

INPUTAGEWEIGHTOXYRUNTIMERSTPULSERUNPULSEMAXPULSE@@;PROC

REG

DATA=D10M04OUTEST=EST;MODELOXY=AGEWEIGHTRUNTIMERUNPULSEMAXPULSERSTPULSE/SELECTION=STEPWISE;MODELOXY=AGEWEIGHTRUNTIMERUNPULSEMAXPULSERSTPULSE/SELECTION=MAXR;RUN;多元线性回归结果：根据逐步回归法的选择标准，依次进入回归模型的自变量是RUNTIME、AGE、RUNPULSE与MAXPULSE。这四个变量的组合可解释84.3%的OXY的变异。模型的形式是吸氧效率OXY的预测值=100.07909519-0.21265570AGE-2.76824065RUNTIME-0.33956528RUNPULSE+0.25535199MAXPULSE此模型的F检验值高达34.90，P<0.0001，是一个显著的回归方程。多元线性回归当某一自变量与其它自变量之间有高度线性相关时，参数的估计值将会不稳定，而且会含偏高的标准误，这个现象称为共线性（Collinearity）或多变量共线性（Multicollinearity）。针对这一个问题，我们可采用COLLIN选项来诊断到底哪些自变量之间有共线性。这个诊断的理论基础来自Belsley、Kuh及Welsch于1980年所发表的论文。多元线性回归诊断的步骤如下第一，将（X’X）矩阵标准化，使其对角线上的值都成为1。若选用COLLINOINT选项，则不包括截距项。第二，计算出（X’X）矩阵的特性根与特性向量。第三，以最大的特性根为分子，其它特性根分别为分母，形成几个不同的比例，这些比例的平方根便是共线性指标。若指标的值较大时，则表示变量之间的共线性情形可能极为严重，在这种情况下，参数的估计值较不准确。多元线性回归示例10M06：诊断示例10M05资料自变量间的相关性。/*SasProg10M06.Sas*/PROC

REG

DATA=d10m05;MODELOXY=RUNTIMEAGEWEIGHTRUNPULSEMAXPULSERSTPULSE/TOL

VIF

COLLIN;RUN;Model语句选项Tol，给出参数估计值的容忍度；选项Vif，将方差膨胀因子随参数估计值一同给出。方差膨胀因子等于容忍度的倒数。多元线性回归CollinearityDiagnostics

NumberEigenvalueConditionIndex

---------------ProportionofVariation---------------InterceptRUNTIMEAGEWEIGHTRUNPULSEMAXPULSERSTPULSE16.949161.000000.000023290.000215130.000163670.000197520.000008580.000006260.0002781120.0192219.015940.001880.021930.174960.005230.000014853.981161E-80.3515630.0151121.448410.000802340.131820.137150.242480.001170.001250.0497740.0091627.548750.005930.631460.030200.168470.001400.001150.2075450.0061433.634350.001770.114480.105760.462710.014650.008160.3647160.0010481.807500.785270.085840.477590.098680.070260.005300.0195270.00017734197.952060.204320.014260.074180.022240.912490.984120.00662共线性比较严重的变量是RUNPULSE与MAXPULSE，此外RUNTIME与RSTPULSE以及WEIGHT与RSTPULSE之间也有中等程度的共线性。岭回归当|X’X|≈0时，我们设想给X’X加上一个正常数矩阵kI（k>0），那么(X’X+

kI)-1接近奇异的可能性就会比(X’X)-1接近奇异的可能性小得多；因此用作为β的估计值应比最小二乘估计稳定，我们称为β的岭回归估计，k为岭参数，常常用较直观的岭迹法确定k值。岭回归当自变量之间存在复共线性时，建立的模型会极不稳定。岭回归模型是通过在自变量自相关矩阵中引入一个很小的岭参数k(0<k<1),并将它加到主对角线元素上,从而降低参数的最小二乘估计复共线特征向量的影响,以保证参数估计更接近真实情况。岭回归示例10L01（郭秀花，赵连伟.SAS6.11版岭回归分析程序设计及其实例分析.数理统计与管理.2001,20(1):41-44,64）：实例选自胡良平主编的“现代统计学与SAS应用”一书。有人在某地抽样调查了29例儿童的血红蛋白与4种微量元素的含量,资料如下,试问可否用4种微量元素(单位都是μmo1/L)钙(X1)、镁(X2)、铁(X3)、铜(X4)来较好地预测血红蛋白(Y,g/L)的含量?数据文件SasData10L01.txt岭回归/*SasProg10L01.Sas*/Datadata01;

infile

'e:\sas\sasdata\sasdata10L01.txt';

Inputyx1-x4;Proc

standard

m=0

s=1

out=data02;run;Proc

reg

data=data02outest=rrr;/*实现岭回归*/

Modely=x1-x4/ridge=0.0to1.0by0.1

collinoint;

Plot

rsq.*np./ridgeplot;/*绘制岭迹图*/Run;Proc

print

data=rrr;Run;STANDARD过程把x1,x2,…,xm,y作标准化转换,转换后的标准化变量均数为0,标准差为1。岭回归结果共线性分析表明，自变量之间的共线性主要表现在X2与X3、X1与X4之间。本例适合采用岭回归进行分析。在不同的岭参数下,各自变量的回归系数保存在数据集rrr中。从岭迹图或从数据集rrr中可以看到，当岭参数大于0.1时，自变量的回归系数趋于平稳(以水平直线为渐进线)，故选岭参数为0.2时的标准化回归方程为：y=6.0663E-16-0.20612x1+0.17949x2+0.65829x3+0.08377x4岭回归结果从岭回归结果（ParameterEstimates）中可以得知：在4种微量元素中，影响血红蛋白的含量达到显著性水平的有：铁及钙。从回归系数的符号可知：铁的含量提高，有助于血红蛋白含量的提高，而钙的吸收量加大后，反面会使血红蛋白含量减少。Nlin过程非线性回归模型（Nonlinearregression）的估计方法有多种思路，常见的有先利用变量变换将模型转换为线性模型，然后进行估计；或者使用泰勒级数展开式进行逐次的线性近似估计，直至达到迭代收敛标准。严格地讲，能够通过变量变换被转化为线性回归方程，并通过对该线性模型的估计以求得参数值的模型只能被称为曲线回归模型（固有线性模型）；直接进行估计的才是非线性回归模型（固有非线性模型）。Nlin过程SAS中主要有两个过程可以实现非线性回归，分别是Nlin过程和Nlmixed过程，后者是非线性混合效应模型。Nlin过程使用最小二乘法或加权最小二乘法实现非线性回归模型的参数估计。与Reg过程相比，Nlin过程中除简单列出参与回归模型的变量外，用户还必须给出方程表达式，定义未知参数名称，并给出所有参数的初始值。由于非线性模型的复杂性，Nlin过程并不保证所有模型都能被拟合出来。Nlin过程Nlin过程的一般格式：PROCNLIN<options>;MODELdependent=expression;PARAMETERS

parameter=values<,...,parameter=values>;BOUNDSinequality<,...,inequality>;BY

variables;DER.parameter=expression;DER.parameter.parameter=expression;ID

variables;OUTPUTOUT=SAS-data-setkeyword=names<,...,keyword=names>;CONTROL

variable<=values><...variable<=values>>;Nlin过程Nlin过程在进行非线性回归模型参数估计时，提供了5种不同的迭代（搜索）方法：Gauss（默认方法），高斯法，又称改良高斯-牛顿法（ModifiedGauss-Newtonmethod）；Marquardt，麦夸特法（Marquardt

method）；Newton，牛顿法（Newtonmethod）。Gradient，梯度法（Gradientmethod），又称最速下降法（Steepestdescentmethod）；Dud，又称错位法（Falseposition）或多元割线法（Multivariatesecant）；可以在选项中使用“method=”指明迭代方法，前三种方法更为稳键。Nlin过程Nlin过程中的model语句：Model语句用来指定所要拟合的非线性回归模型。格式为：Modeldependent=expression;此处的Model语句必须给出模型完整的表达式，比如：Modely=exp(a+b*x);指数模型Modely=a+b*cos(2*pi*t/12+c);余弦模型Modely=b0*(1-exp(-b1*x);负指数增长曲线未知参数可以自行命名，表达式中可以使用SAS函数和常数。Nlin过程Nlin过程中的Parameters语句：Parameters语句也可以简写作Parms。由于Model语句中出现了未知参数，因此需要Parameters语句给定参数初始值。通常，每个参数只需设定一个初始值即可。如：parmsa=1b=1;给定离散的多个值：m1,m2,…,mn；搜索整数区间：mTon；从m到n，间隔I的等差序列：mTonByi；离散数值和整数区间：m1,m2Tom3。Nlin过程Nlin过程中可以使用类似于Data

Step中的相应赋值语句，经实现较为复杂的非线性模型，如下面的程序：Procnlin;parmsb0=0b1=0.01;temp=exp(b1*x);modely=b0*(1-temp);显然，上述模型等价于y=b0*(1-exp(b1*x))，但表达式更为简洁，并且拟合的速度更快。除赋值语句外，绝大多数Data

Step中可用的编程语句，如判断、循环、goto语句等都可在这里使用。非线性回归分析实例10N01：收集我国1940-1998年间不同月份的男婴死亡率（‰），试为全年各月份的死亡率波动拟合适当的模型。月份123456死亡率63.855.854.557.058.553.5月份789101112死亡率57.367.965.468.567.462.2非线性回归分析首先，建立数据集，研究死亡率与月份之间的关系（散点图）。/*SasProg10n01.Sas*/dataD10n01;

inputxy@@;cards;163.8255.8354.5457.0558.5653.5757.3867.9965.41068.51167.41262.2;proc

gplot;

ploty*x='+';run;非线性回归分析非线性回归分析由图可见，死亡率3-6月间达到最低，而8-11月份达到最高，在全年范围内呈现一种有规律的波动趋势。试着拟合如下模型：y=b0+b1*cos(2*3.1415926/12*x+b2)其中：b0,b1,b2为未知数；

2*3.1415926/12*x用于将月份转变为弧度。非线性回归分析/*SasProg10n02.Sas*/proc

nlin

data=d10n01;

modely=b0+b1*cos(2*3.1415926/12*x+b2);

parametersb0=60b1=6b2=4;run;选项maxiter=设定迭代次数，默认次数100；选项Converge=设定Bates和Watts提出的收敛测量指标，默认收敛标准10-5。选项ConvergeObj=设定以损失函数的改变量作为收敛标准，默认值10-8。非线性回归分析TheNLINProcedureDependentVariableyMethod:Gauss-NewtonIterativePhaseSumofIterb0b1b2Squares060.00006.00004.00001028.5160.9833-6.47223.866296.4181260.9833-6.30724.122680.0095360.9833-6.52184.123579.7329460.9833-6.52184.123479.7329依变量为y迭代方法：高斯-牛顿法每次迭代过程中的参数估计值以及残差平方和的变化。可以看出，随着迭代次数的增加，各参数估计值越来越稳定，残差平方和越来越小。输出的第一部分：非线性回归分析模型的第二部分给出收敛的状态（Convergencecriterionmet.合适的）。输出的第三部分给出模型拟合情况汇总，迭代方法、迭代次数、相应的模型收敛标准等。（可忽略）输出的第四部分是近似方差分析结果。相应的原假设是：所拟合的模型对依变量的预测无贡献。显示，p=0.0016<0.05，可以认为模型中至少有一些自变量对依变量是有作用的。输出的第五部分是模型中未知参数的点估计和区间估计。输出的第六部分为各参数的相关系数阵。本例三个未知参数间不存在相关关系。非线性回归分析实例10N03（张录达编.数据拟合.中国农业大学自编教材）：植物病毒对植物叶片的伤害（一般取枯斑数）与植物病毒的浓度之间是一种非线性关系。一般在低的、中等病毒浓度下二者之间的关系可用下式表示：其中，y为平均侵染枯斑数，v为病毒浓度。N、C为未知数。病毒浓度较高时，植物病毒侵-稀度关系的数学模型为：其中，y为半叶侵染的的枯斑数，v为病毒浓度。非线性回归分析实例10N03（张录达编.数据拟合.中国农业大学自编教材）：不同浓度下心叶烟接种烟草花叶病毒（TMV），在半叶所得平均枯斑数据。试以的极小化为目标，采用Marquardt法进行植物侵染-稀度模型的数据拟合。非线性回归分析实例10N03（张录达编.数据拟合.中国农业大学自编教材）：病毒浓度(mg/ml)组别12345671.0E-14652341492541963161923.5E-22921851112161642522211.0E-222811970.615691.51901333.5E-316083.070.010881.613988.21.0E-376.549.830.059.637.572.254.83.5E-459.427.510.933.626.433.327.31.0E-420.08.185.514.513.117.510.23.5E-58.43.253.858.067.089.092.941.0E-5.2.441.314.261.644.23.非线性回归分析/*SasProg10n03.Sas*/datad10n03;

inputxy1-y7;av=mean(ofy1-y7);cards;1.0E-14652341492541963161923.5E-22921851112161642522211.0E-222811970.615691.51901333.5E-316083.070.010881.613988.21.0E-376.549.830.059.637.572.254.83.5E-459.427.510.933.626.433.327.31.0E-420.08.185.514.513.117.510.23.5E-58.43.253.858.067.089.092.941.0E-5.2.441.314.261.644.23.;Mean是SAS提供的求样本平均数的函数。非线性回归分析proc

nlin

method=marquardt;

parmsb1=30b2=4000;modely1=b1*log(1+b2*x);proc

nlin

method=marquardt;

parmsb1=10b2=2000;modely2=b1*log(1+b2*x);proc

nlin

method=marquardt;

parmsb1=10b2=5000;modely3=b1*log(1+b2*x);proc

nlin

method=marquardt;

parmsb1=10b2=6000;modely4=b1*log(1+b2*x);proc

nlin

method=marquardt;

parmsb1=60b2=8000;modely5=b1*log(1+b2*x);proc

nlin

method=marquardt;

parmsb1=10b2=5000;modely6=b1*log(1+b2*x);proc

nlin

method=marquardt;

parmsb1=20b2=1000;modely7=b1*log(1+b2*x);proc

nlin

method=marquardt;

parmsb1=20b2=1000;modelav=b1*log(1+b2*x);run;非线性回归分析组别方程显著性（ApproxPr>F）1y1=87.6185ln(1+1342.2v)<.00012y2=45.2684ln(1+1618.2v)<.00013y3=26.8680ln(1+2068.6v)<.00014y4=45.1564ln(1+2998.1v)<.00015y5=37.2329ln(1+1857.2v)<.00016y6=55.5436ln(1+2869.0v)<.00017y7=37.4928ln(1+3499.1v)<.0001Meanav=47.2415ln(1+2204.4v)<.0001对数曲线拟合实例10D01（阮桂海等编著.SAS统计分析实用大全.北京:清华大学出版社.2003,6.第1版:408-413）：一组产量（y）与劳动力（L）、资金（K）的数据，试拟合非线性经济模型：参数：B0：截距；D：分布参数；C：斜率，即效率参数；r：替代参数。数据文件：SasData10d01.txt对数曲线拟合/*SASProg10d01.Sas*/title

'非线性回归分析-对数曲线回归';datad10d01;

infile

'e:\sas\sasdata\sasdata10d01.txt';

inputkly;proc

nlin

best=100

method=dud;

modely=b0+c*log(d*l**r+(1-d)*k**r);

parmsb0=1c=-1d=0.5r=-1;run;选项best=n，只列出整个搜索过程前n个最小的残差平方和。对数曲线拟合输出结果的解释：第一部分、第二部分：模型的初始参数系数值。第三部分：迭代了37步，最小平方和已经不能再下降了（达到收敛标准）。第四部分：

NOTE:Convergencecriterionmet.收敛标准是适合的第五部分：模型拟合情况汇总。第六部分：方差分析Regressionsumofsquares(RegressionSS)：已被解释的回归平方和。Residualsumofsquares(ResidualSS)：未被解释的回归平方和。Uncorrectedtotalsumofsquares：依变量的总平方和。Correctedtotalsumofsquares：偏离均值的平方和。Rsquare=1-ResidualSS/CorrectedSS=1-1.6492/58.2089=0.9717，意味着：能被模型解释的方差点总方差的97.17%。Logistic过程如果依变量为分类变量，则不符合一般回归分析模型的要求，可用Logistic回归来分析。Logistic回归按反应变量的类型分为：两分类的Logistic回归多分类有序反应变量的Logistic回归多分类无序反应变量的Logistic回归按照设计类型可分为：非条件Logistic回归，即研究对象未经过配对条件Logistic回归，即研究对象为1：1或1：m配对简单的Logistic回归需调用SAS中LOGISTIC过程完成，一些较复杂的则需要调用CATMOD过程来实现。Logistic过程Logisticregression回归中，可直接计算一个事件的发生概率：有多个自变量的Logisticregression模型，可以写成：式中，Z中线性结合模型，即事件没有发生的概率可以写成：

Prob(noevent)=1-Prob(event)Logistic过程优势（odds）：发生概率（p）与不发生概率（1-p）之比。对odds取自然对数：Logistic过程Logistic过程语法：PROCLOGISTIC<options>;BYvariables;CLASSvariable<(v-options)><variable<(v-options)>...></v-options>;CONTRAST'label'effectvalues<,...effectvalues></options>;FREQvariable;MODELresponse=<effects></options>;MODELevents/trials=<effects></options>;OUTPUT<OUT=SAS-data-set><keyword=name...keyword=name>/<option>;<label:>TESTequation1<,...,<equationk>></option>;UNITSindependent1=list1<...independentk=listk></option>;WEIGHTvariable</option>;Logistic过程Logistic过程用最大似然法对应变量拟合一个Logistic模型。除了ProcLogistic和Model语句为必需，其他都可省略。Logistic过程示例10C01（GeoffDer,BrianS.Everitt.AHandbookofStatisticalAnalysesusingSAS(SECONDEDITION).CHAPMAN&HALL/CRCPress.2002）：精神病调查问卷，内容包括得分、性别、病例数、非病例数。数据文件SasData10c01.txt首先，建立数据集ghq（GeneralHealthQuestionnaire）；生成总调查数(Total)和发病率(Prcases)变量；绘制得分与发病率散点图。Logistic过程/*SasProg10c01.Sas*/dataghq;

infile

'e:\sas\sasdata\sasdata10c01.txt';

inputghqsex$casesnoncases;total=cases+noncases;prcase=cases/total;proc

gplot

data=ghq;

plotprcase*ghq;run;结果显示，随着得分的增大，患病的可能性增大。Logistic过程Logistic过程/*SasProg10c02.Sas*/proc

reg

data=ghq;modelprcase=ghq;output

out=routp=rpred;run;proc

logistic

data=ghq;modelcases/total=ghq;output

out=loutp=lpred;run;我们首先实施线性拟合，然后再用Logistic模型拟合；分别将结果输出到rout、lout数据集，数据集中包含预测值（rpred和lpred变量，患病可能性）。Logistic过程有两种形式的Model语句，本例使用了events/trials语法，两个变量被一个斜杠分开，指定一个“二择一”事件。Logistic过程模型统计量（ModelFitStatistics）给出了仅由截距拟合的模型似然函数和基于截距和协变量共同拟合的模型似然函数的三种判断标准：AIC、SC、-2LogL（即Likelihood）。无效假设为β=0的假设测验结果（TestingGlobalNullHypothesis:BETA=0）：卡方检验概率都小于0.01，所以模型具有极显著意义。最大似然估计分析（AnalysisofMaximumLikelihoodEstimates）：得到对数线性判别函数logit(p)=-2.7107+0.7360ghqLogistic过程优势率估计（OddsRatioEstimates）：分为点估计和95%置信区间估计两部分。本例可以看出，ghq优势增长一个单位患病的机率增长2.088（1.734-2.513）倍。最后输出的内容为预测概率与观测依变量之间的关联性（AssociationofPredictedProbabilitiesandObservedResponses）。Logistic过程/*SasProg10c03.Sas*/datalrout;setrout;setlout;proc

sort

data=lrout;byghq;proc

gplot

data=lrout;plot(rpredlpredprcase)*ghq/overlay;run;原始数据及两种回归分析方法预测值数据集合并。然后，共同绘散点图。Logistic过程Logistic过程从两种回归预测值集合，我们就能看出，线性回归模型是不适合的。原因是响应变量是一种可能性（概率），取值范围（0，1），线性模型的预测值超出了这个范围，并且拟合也不好。Logistic模型的预测值都落在0-1之间，对观测数据拟合也较好。Logistic过程/*SasProg10c04.Sas*/proc

logistic

data=ghq;classsex;modelcases/total=sexghq;run;对本例进行扩展，将“性别”因素考虑在内进行分析。Class语句指定分类变量；Model语句指定效应模型，与glm过程相似。Logistic过程结果显