r - corr.test arguments imply differing number of rows - Stack Overflow
Join Stack Overflow to learn, share knowledge, and build your career.
or sign in with
I have seen this error multiple times in different projects and I was wondering if there is a way to tell which line caused the error in general?
My specific case:
http://archive.ics.uci.edu/ml/machine-learning-databases/00275/
#using the bike.csv
data&-read.csv("PATH_HERE\\Bike-Sharing-Dataset\\day.csv",header=TRUE)
require(psych)
corr.test(data)
data&-data[,c("atemp","casual","cnt","holiday","hum","mnth","registered",
"season","temp","weathersit","weekday","windspeed","workingday","yr")]
data[data=='']&-NA
#View(data)
require(psych)
cors&-corr.test(data)
returns the error:
Error in data.frame(lower = lower, r = r[lower.tri(r)], upper = upper,
arguments imply differing number of rows: 0, 91
2,28552154
It works for me
& #using the bike.csv
& data &- read.csv("day.csv",header=TRUE)
& require(psych)
& corr.test(data)
Error in cor(x, use = use, method = method) : 'x' must be numeric
& data &- data[,c("atemp","casual","cnt","holiday","hum","mnth","registered",
"season","temp","weathersit","weekday","windspeed","workingday","yr")]
& data[data==''] &- NA
& #View(data)
& require(psych)
& cors &- corr.test(data)
Call:corr.test(x = data)
Correlation matrix
atemp casual
cnt holiday
mnth registered season
0.05 -0.08
-0.07 -0.10
0.05 -0.07
1.00 -0.02
-0.01 -0.03
-0.08 -0.10
registered
-0.11 -0.09
weathersit -0.12
-0.25 -0.30
0.02 -0.12
-0.10 -0.05
-0.17 -0.23
0.01 -0.25 -0.21
-0.23 -0.16
workingday
0.02 -0.01
0.01 -0.11
weathersit weekday windspeed workingday
0.02 -0.11
registered
weathersit
0.06 -0.05
0.04 -0.01
-0.02 -0.01
workingday
Sample Size
Probability values (Entries above the diagonal are adjusted for multiple tests.)
atemp casual
cnt holiday
hum mnth registered season temp
1.00 0.01 0.00
1.00 1.00 0.04
1.00 0.28 0.00
0.00 1.00 1.00
0.67 0.00 0.00
0.60 0.00 0.00
registered
0.00 0.01 0.00
0.78 0.00 0.00
0.44 0.00 0.00
weathersit
0.35 0.00 0.24
0.01 0.16 0.80
0.87 0.00 0.00
workingday
0.00 0.51 0.87
0.83 0.00 0.96
weathersit weekday windspeed workingday
registered
weathersit
workingday
To see confidence intervals of the correlations, print with the short=FALSE option
4,73211431
It works for me:::
rm(list=ls())
# http://archive.ics.uci.edu/ml/machine-learning-databases/00275/
#using the bike.csv
day &- read.csv("Bike-Sharing-Dataset//day.csv")
require(psych)
day&-day[,c("atemp","casual","cnt","holiday","hum","mnth","registered",
"season","temp","weathersit","weekday","windspeed","workingday","yr")]
day[day=='']&-NA
require(psych)
corr.test(day)
# corr.test(day)
# Call:corr.test(x = day)
# Correlation matrix
# atemp casual
cnt holiday
mnth registered season
temp weathersit weekday windspeed workingday
0.05 -0.08
-0.07 -0.10
0.05 -0.07
1.00 -0.02
-0.01 -0.03
-0.08 -0.10
0.02 -0.11
# registered
-0.11 -0.09
# weathersit -0.12
-0.25 -0.30
0.02 -0.12
0.06 -0.05
-0.10 -0.05
0.04 -0.01
# windspeed
-0.17 -0.23
0.01 -0.25 -0.21
-0.23 -0.16
-0.02 -0.01
# workingday
0.02 -0.01
0.01 -0.11
# Sample Size
# Probability values (Entries above the diagonal are adjusted for multiple tests.)
# atemp casual
cnt holiday
hum mnth registered season temp weathersit weekday windspeed workingday
1.00 0.01 0.00
1.00 1.00 0.04
1.00 0.28 0.00
0.00 1.00 1.00
0.67 0.00 0.00
0.60 0.00 0.00
# registered
0.00 0.01 0.00
0.78 0.00 0.00
0.44 0.00 0.00
# weathersit
0.35 0.00 0.24
0.01 0.16 0.80
# windspeed
0.87 0.00 0.00
# workingday
0.00 0.51 0.87
0.83 0.00 0.96
# To see confidence intervals of the correlations, print with the short=FALSE option
Your Answer
Sign up or
Sign up using Google
Sign up using Facebook
Post as a guest
Post as a guest
By posting your answer, you agree to the
Not the answer you're looking for?
Browse other questions tagged
Upcoming Events
ends Mar 27
Stack Overflow works best with JavaScript enabledR语言的各种检验_百度文库
两大类热门资源免费畅读
续费一年阅读会员，立省24元！
R语言的各种检验
你可能喜欢r - Can corr.test be used for any dataframe? - Stack Overflow
Join Stack Overflow to learn, share knowledge, and build your career.
or sign in with
I am trying to find get correlations and p-values between variables in a dataframe (df1) using corr.test in the psych package.
The variables in the
dataframe are all integers and there is no NAs.
But when I run the corr.test(df1), there is always a error message.
Error in data.frame(lower = lower, r = r[lower.tri(r)], upper = upper,
arguments imply differing number of rows: 0, 28
I tried to run the example (corr.test(sat.act)) in the psych package and there is no error.
I am new to R, can someone tell me what is wrong with the dataframe.
& head(df1)
S1.pre S2.pre S1.post S2.post V1.pre V2.pre V1.post V2.post
& dput(head(df1))
structure(list(S1.pre = c(21L, 15L, 18L, 17L, 15L, 17L), S2.pre = c(31L,
26L, 27L, 31L, 26L, 28L), S1.post = c(25L, 21L, 23L, 18L, 16L,
16L), S2.post = c(35L, 29L, 31L, 39L, 29L, 28L), V1.pre = c(7L,
13L, 8L, 13L, 26L, 2L), V2.pre = c(1L, 11L, 2L, 11L, 15L, 4L),
V1.post = c(19L, 16L, 3L, 15L, 32L, 2L), V2.post = c(4L,
14L, 3L, 14L, 20L, 7L)), .Names = c("S1.pre", "S2.pre", "S1.post",
"S2.post", "V1.pre", "V2.pre", "V1.post", "V2.post"), row.names = c(NA,
6L), class = "data.frame")
& sapply(df1, class)
"integer" "integer" "integer" "integer" "integer" "integer" "integer" "integer"
I contacted William Revelle - author of the psych package and here is what he said:
Unfotunately you found a bug introduced into 1.4.3.
1.4.4 will go out to Cran this weekend.
In the meantime you can get the fix at
(choose source from other repository if you are using a mac) or
and get the zip file if you are using a PC.
Otherwise, wait until next week.
Sorry about the problem.
It will still work as long as you have unequal number of subjects or some missing data.
Your Answer
Sign up or
Sign up using Google
Sign up using Facebook
Post as a guest
Post as a guest
By posting your answer, you agree to the
Not the answer you're looking for?
Browse other questions tagged
Upcoming Events
ends Mar 27
Stack Overflow works best with JavaScript enabled在 SegmentFault，学习技能、解决问题
每个月，我们帮助 1000 万的开发者解决各种各样的技术问题。并助力他们在技术能力、职业生涯、影响力上获得提升。
问题对人有帮助，内容完整，我也想知道答案
问题没有实际价值，缺少关键内容，没有改进余地
做相似度计算的时候经常会用到皮尔逊相关系数（Pearson Correlation Coefficient），那么应该如何理解该系数？其数学本质、含义是什么？
相关公式：
相关链接：
答案对人有帮助，有参考价值
答案没帮助，是错误的答案，答非所问
回报一下 @xhinking皮尔逊相关系数理解有两个角度其一, 按照高中数学水平来理解, 它很简单, 可以看做将两组数据首先做Z分数处理之后, 然后两组数据的乘积和除以样本数Z分数一般代表正态分布中, 数据偏离中心点的距离.等于变量减掉平均数再除以标准差.(就是高考的标准分类似的处理)标准差则等于变量减掉平均数的平方和,再除以样本数,最后再开方.所以, 根据这个最朴素的理解,我们可以将公式依次精简为:其二, 按照大学的线性数学水平来理解, 它比较复杂一点,可以看做是两组数据的向量夹角的余弦.皮尔逊相关的约束条件从以上解释, 也可以理解皮尔逊相关的约束条件:1 两个变量间有线性关系2 变量是连续变量3 变量均符合正态分布,且二元分布也符合正态分布4 两变量独立在实践统计中,一般只输出两个系数,一个是相关系数,也就是计算出来的相关系数大小,在-1到1之间;另一个是独立样本检验系数,用来检验样本一致性.先举个手算的例子使用维基中的例子:例如，假设五个国家的国民生产总值分别是1、2、3、5、8（单位10亿美元），又假设这五个国家的贫困比例分别是11%、12%、13%、15%、18%。创建2个向量.(R语言)x&-c(1,2,3,5,8)y&-c(0.11,0.12,0.13,0.15,0.18)按照维基的例子,应计算出相关系数为1出来.我们看看如何一步一步计算出来的.x的平均数是:3.8
y的平均数是0.138
所以,sum((x-mean(x))*(y-mean(y)))=0.308用大白话来写就是:(1-3.8)*(0.11-0.138)=0.0784
(2-3.8)*(0.12-0.138)=0.0324
(3-3.8)*(0.13-0.138)=0.0064
(5-3.8)*(0.15-0.138)=0.0144
(8-3.8)*(0.18-0.138)=0.17640.4+0.4+0.同理, 分号下面的,分别是:sum((x-mean(x))^2)=30.8
sum((y-mean(y))^2)= 0.00308用大白话来写,分别是:(1-3.8)^2=7.84
(2-3.8)^2=3.24
(3-3.8)^2=0.64
(5-3.8)^2=1.44
(8-3.8)^2=17.64
#平方7.84+3.24+0.64+1.44+17.64=30.8同理,求得:sum((y-mean(y))^2)= 0.00308然后再开平方根,分别是:30.8^0.5=5.549775
0.=0.用分子除以分母,就计算出最终结果:0.308/(5..再举个简单的R语言例子(R在这里下载: )假设有100人, 一组数据是年龄,平均年龄是35岁,标准差是5岁;另一组数据是发帖数量,平均帖子数量是45份post,标准差是8份帖子.假设这两组都是正态分布.我们来求这两者的皮尔逊相关系数,R脚本如下:& x&-rnorm(n=100,mean=35,sd=5)
#创建一组平均数为35,标准差为5,样本数为100的随机数
& y&-rnorm(n=100,mean=45,sd=8) #创建一组平均数为45,标准差为8,样本数为100的随机数
cor.test(x,y,method=&pearson&) #计算这两组数的相关,并进行T检验然后R输出结果为:Pearson's product-moment correlation
t = -0.0269, df = 98, p-value = 0.9786
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
-0..1938019
sample estimates:
-0.当然,这里是随机数.也可以用非随机的验证一下计算.皮尔逊相关系数用于网站开发直接将R与Ruby关联起来调用很简单,仿照上述例子:cor(x,y)就输出系数结果了.有这么几个库可以参考:说明, 以上为ruby调用库. pythone程序员可以参考: Rpy ()简单的相关系数的分类0.8-1.0 极强相关0.6-0.8 强相关0.4-0.6 中等程度相关0.2-0.4 弱相关0.0-0.2
极弱相关或无相关ps : 这个网站开发者不要再次发明轮子,本来用markdown语法写作很爽,结果又不得不花时间来改动.请考虑尽快支持Markdown语法.
答案对人有帮助，有参考价值
答案没帮助，是错误的答案，答非所问
皮尔森相关系数的就是 $$\frac{x和y的协方差}{(x的标准差 * y的标准差)}$$
判断两组数的线性关系程度。
答案对人有帮助，有参考价值
答案没帮助，是错误的答案，答非所问
玩一玩这个就知道了
同步到新浪微博
分享到微博？
关闭理由：
删除理由：
忽略理由：
推广（招聘、广告、SEO 等）方面的内容
与已有问题重复（请编辑该提问指向已有相同问题）
答非所问，不符合答题要求
宜作评论而非答案
带有人身攻击、辱骂、仇恨等违反条款的内容
无法获得确切结果的问题
非开发直接相关的问题
非技术提问的讨论型问题
其他原因（请补充说明）
我要该，理由是：
在 SegmentFault，学习技能、解决问题
每个月，我们帮助 1000 万的开发者解决各种各样的技术问题。并助力他们在技术能力、职业生涯、影响力上获得提升。SAS的基本统计功能
&&您的位置：
&&&&SAS的基本统计功能
SAS的基本统计功能
前面我们已经看到了SAS的编程计算、数据管理能力、数据汇总、数据探索分析能力。这一章我们讲如何用SAS进行基本的统计检验、线性回归、方差分析、列联表检验等基本统计分析。我们既使用SAS语言编程，也使用SAS/INSIGHT的菜单界面。
对单个变量，我们可能需要作正态性检验、两独立样本均值相等的检验、成对样本均值相等的检验。
在PROC UNIVARIATE语句中加上NORMAL选项可以进行正态性检验。例如，我们要检验SASUSER.GPA 中GPA是否服从正态分布，只要用如下UNIVARIATE过程：
proc univariate data=sasuser.
结果（部分）如下：
Univariate Procedure
Variable=GPA
College Grade Point Average
其中W:Normal为Shapiro-Wilk正态性检验统计量，Pr&W为检验的显著性概率值（p值）。当N≤2000时正态性检验用Shapiro-Wilk统计量，N&2000时用Kolmogorov
D统计量。我们可以看到，p值很小，所以在0.05水平（或0.10水平）下应拒绝零假设，即认为GPA分布非正态。
在SAS/INSIGHT中为了检验GPA的分布，先选“Analyze | Distribution”菜单打开GPA 变量的分布窗口，然后选“Curves
| Test for Distribution”菜单。除了可以检验是否正态分布外还可以检验是否对数正态、指数分布、Weibull分布。
假设我们有两组样本分别来自两个独立总体，需要检验两个总体的均值或中心位置是否一样。如果两个总体都分别服从正态分布，而且方差相等，可以使用两样本t检验过程TTEST
比如，我们要检验SASUSER.GPA数据集中男生和女生的SATM分数是否具有相等的平均值，只要用如下程序：
proc ttest data=sasuser.
过程中用CLASS语句指定分组变量，用VAR语句指定要比较的变量。结果如下：
TTEST PROCEDURE
Variable: SATM
Math SAT Score
-----------------------------------------------------------------------------
---------------------------------------
For H0: Variances are equal, F' = 1.03
DF = (144,78)
Prob&F' = 0.9114
结果有三个部分：两个总体的SATM简单统计量，两样本均值的检验，以及两样本方差是否相等的检验。标准的两样本t检验要求两总体方差相等，所以第三部分结果检验两样本方差是否相等。如果检验的结果为相等，则可使用精确的两样本t检验，看第二部分结果的Equal那一行。如果方差检验的结果为不等，则只能使用近似的两样本t检验，看第二部分结果的Unequal
那一行。这里我们看到方差检验的p值为0.9114不显著，所以可以认为方差相等，所以我们看Equal 行，p值为0.水平下是显著的，所以应认为男、女生的SATM分数有显著差异，女生分数要高。
上面的检验中对立假设是两组的均值不等，所以检验是双边的，p值的计算公式为Pr(t 分布随机变量绝对值&计算得到的t统计量的绝对值)。如果要进行单边的检验，比如对立假设为女生分数高于男生分数（右边），则p值为Pr(t分布随机变量&计算得到的t统计量)
，当计算得到的t统计量值为正数时（现在t=4.0）此单边p值为双边p值的一半，当计算得到的t统计量为负数时肯定不能否定零假设。检验左边时恰好相反。
如果我们希望检验男、女生的GPA分数则无法使用两样本t检验，因为检验女生的GPA 样本的正态性发现它非正态。这种情况下我们可以使用非参数检验。检验两独立样本的位置是否相同的非参数检验有Wilcoxon秩和检验。我们用NPAR1WAY过程加Wilcoxon选项可以进行这种检验。见下例：
proc npar1way data=sasuser.
其CLASS语句和VAR与TTEST过程相同。结果如下：
N P A R 1 W A Y
P R O C E D U R E
Wilcoxon Scores (Rank Sums) for Variable GPA
Classified by Variable SEX
463.429146
110.810345
463.429146
115.601266
Average Scores Were Used for Ties
Wilcoxon 2-Sample Test (Normal Approximation)
(with Continuity Correction of .5)
Z = 0.527589
Prob & |Z| = 0.5978
T-Test Approx. Significance = 0.5983
Kruskal-Wallis Test (Chi-Square Approximation)
CHISQ = 0.27949
Prob & CHISQ = 0.5970
结果分为四部分：两样本的秩和的有关统计量，Wilcoxon两样本检验的结果，t检验的近似显著性，Kruskal-wallis检验结果。我们只要看Wilcoxon检验的p值Prob
& |Z| = 0.5978 ，检验结果不显著，可认为男、女生的GPA分数在0.05水平下无显著差异。
SAS/INSIGHT中未提供两独立样本检验的功能。
我们在现实中经常遇到两个总体是相关的测量结果的比较，比如，考察同一组人在参加一年的长跑锻炼前后的心率有无显著差异。这时，每个人一年前的心率和一年后的心率是相关的，心率本来较快的人锻炼后仍相对于其它人较快。所以，检验这样的成对总体的均值不能使用两样本t检验的方法，因为独立性条件不再满足。这时，我们可以检验两个变量间的差值的均值是否为零，这等价于检验两组测量值的平均水平有无显著差异。
检验单个样本的均值是否为零只要使用UNIVARIATE过程，在UNIVARIATE过程的矩部分给出了均值为零的t检验和符号检验、符号秩检验的结果。例如，我们想知道SATM和SATV这两门考试的成绩有无显著差异（SATM平均值为595.3，SASTV平均值为504.6，我们希望知道差异是否显著）。因为这两个成绩是同一个学生的成绩，所以它们之间是相关的（学得好的学生两科一般都好，学得差的一般两科都差），不能用独立两样本的t检验，但可以计算两变量间的差DMV
＝SATM－SATV，检验差值变量的均值是否为零。如果否定，则可认为SATM和SATV的平均值有显著差异。
为此，我们先用一个数据步计算差值，然后对差值变量用UNIVARIATE过程进行分析就可以得到结果。程序如下：
set sasuser.
dmv = satm -
proc univariate data=
结果（部分）如下：
Univariate Procedure
Variable=DMV
我们只要看其中的三个检验：T: Mean=0是假定差值变量服从正态分布时检验均值为零的t 统计量值，相应的p值Pr&|T|为0.水平下是显著的，所以可认为两科分数有显著差异。M(Sign)是非参数检验符号检验的统计量，其p值Pr&=|M|为0.水平下是显著的，结论不变。Sgn
Rank是非参数检验符号秩检验的统计量，其p值Pr&=|S|为0.0001 在0.05水平下是显著的，结论不变。所以这三个检验的结论都是两科成绩有显著差异。
如果t检验对立假设是单边的，其p值算法与上面讲的两样本t检验p值算法相同。
在SAS/INSIGHT中比较成对样本均值的显著差异，同样是先计算两变量的差值变量（在“Edit | Variables | Other”菜单中，指定两个变量，指定两个变量间的计算为减法，则可以生成差值变量，可以用数据窗口菜单的“Define
Variables”改变量名），然后对此差值变量选“Analyze | Distribution”，选“Tables | Location
Tests”并选中t检验、符号检验和符号秩检验即可在分布窗口显示结果。
本节先讲述如何用SAS/INSIGHT进行曲线拟合，然后进一步讲如何用SAS/INSIGHT进行线性回归，简单介绍SAS/INSIGHT的广义线性模型拟合，最后介绍如何用编程进行回归分析。
两个变量Y和X之间的相关关系经常可以用一个函数来表示，一元函数可以等同于一条曲线，实际工作中经常对两个变量拟合一条曲线来近似它们的相关关系。最基本的“曲线”
是直线，还可以用多项式、样条函数、核估计和局部多项式估计。其模型可表示为
例如，我们要研究SASUSER.CLASS数据集中学生体重与身高之间的相关关系。为此，我们可以先画出两者的散点图（Analyze | Scatter
plot）。从图中可以看出，身高越高的人一般体重越重。我们可以把体重作为因变量、身高作为自变量拟合一条回归直线，只要选 “Analyze |
Fit (Y X)”，并选体重为Y变量，身高为X变量，即可自动拟合出一条回归直线，见图 1。窗口中还给出了拟合的模型方程、参数估计、诊断信息等，我们在下一小节再详细介绍。
在拟合了直线后，为拟合多项式曲线，只要选“Curves | Polynomial”，然后输入阶次（Degree(Polynomial)），就可以在散点图基础上再加入一条多项式曲线。对于本例，我们看到二次多项式得到的曲线与直线差别很小，所以用二次多项式拟合没有优势。还可以试用三次、四次等多项式。为了改变阶次还可以使用拟合窗口中的多项式阶次滑块（Parametric
Regression Fit中的Degree(Polynomial)）。
样条曲线是一种非参数回归的曲线拟合方法。光滑样条为分段的三次多项式，曲线在每一段内是一个三次多项式，在两段的连接点是连续、光滑的。为拟合样条曲线，只要选“Curves
| Spline”，使用缺省的GCV准则（广义交叉核实）来选取光滑系数（光滑系数c越大，得到的曲线越光滑，但拟合同时变差，光滑系数c小的时候得到的曲线较曲折，而拟合较好），就可以在散点图的基础上画出样条曲线。可以用光滑系数c的滑块来调整曲线的光滑程度/拟合优度。对于本例，GCV准则得到的样条曲线与回归直线几乎是重合的，说明直线拟合可以得到满意的结果。
核估计是另一种非参数回归的曲线拟合方法。它定义了一个核函数 ，例如使用标准正态分布密度曲线，然后用如下公式估计经验公式 ：
其中 为光滑系数， 越大得到的曲线越光滑。为了画核估计曲线，只要选“Curves | Kernel”，权重函数使用缺省的正态核，选取光滑系数的方法采用缺省的GCV法，就可以把核估计图附加到散点图上。本例得到的核估计曲线与回归直线、样条曲线有一定差别。可以手动调整光滑系数
的值，可以看到，当 过大时曲线不仅变光滑而且越来越变水平，因为这时的拟合值基本是一个常数，这与样条曲线的情形不同，样条曲线当 增大时曲线变光滑但不趋向于常数（水平线）。
局部多项式估计（Loess）是另一种非参数回归的曲线拟合方法。它在每一自变量值处拟合一个局部多项式，可以是零阶、一阶、二阶，零阶时与核估计相同。SAS/INSIGHT缺省使用一阶（线性）局部多项式。改变Loess的系数alpha可以改变曲线的光滑度。alpha增大时曲线变光滑，而且使用一阶或二阶多项式时曲线不会同时变水平。
固定带宽的局部多项式是另一种局部多项式拟合方法。它有一个光滑系数c。
上面我们已经看到，用菜单“Analyze | Fit (Y X)”就可以拟合一条回归直线，这是对回归方程
的估计结果。这样的线性回归可以推广到一个因变量、多个自变量的情况。线性模型写成矩阵形式为
下面列出了线性模型中常用的一些量和结论：
为 因变量向量
为 矩阵，一般第一列元素全是1，代表截距项
为 未知参数向量
为 随机误差向量，元素独立且方差为相等的 （未知）。
正常情况下，系数的估计为
拟合值（或称预报值）为
空间内向 的列张成的线性空间 投影的投影算子矩阵，叫做“帽子”矩阵。
拟合残差为
残差平方和为
误差项方差的估计为（要求设计阵 满秩）均方误差（MSE）
在线性模型的假设下，若设计阵 满秩，
和 分别是 和 的无偏估计，系数估计的方差阵 。
判断回归结果优劣的一个重要指标为复相关系数平方（决定系数） （其中 ），它代表在因变量的变差中用模型能够解释的部分的比例，所以 越大说明模型越好。
例如，我们在“Fit (Y X)”的选择变量窗口选Y变量（因变量）为体重（WEIGHT），选X 变量（自变量）为身高（HEIGHT）和年龄（AGE），则可以得到体重对身高、年龄的线性回归结果。下面对基本结果进行说明。
回归基本模型：
Response Distribution:
Link Function:
回归模型方程：
Model Equation
拟合概况：
Summary of Fit
Mean of Response
其中Mean of Response为因变量（Response）的均值，Root MSE叫做根均方误差，是均方误差的平方根，R-Square即复相关系数平方，Adj
R-Sq为修正的复相关系数平方，其公式为 ，其中 当有截距项时取1，否则取0，这个公式考虑到了自变量个数 的多少对拟合的影响，原来的 随着自变量个数的增加总会增大，而修正的 则因为 对它有一个单调减的影响所以 增大时修正的 不一定增大，便于不同自变量个数的模型的比较。
方差分析表：
Analysis of Variance
Sum of Squares
Mean Square
这是关于模型是否成立的最重要的检验。它检验的是 ：模型中所有斜率项系数都等于零，这等价于说自变量的线性组合对因变量没有解释作用。它依据的是一个标准的方差分解，把因变量的总离差平方和（C
Total）分解为能用模型解释的部分（Model）与不能被模型解释的部分（随机误差，Error）之和，如果能解释的部分占的比例大就否定 。F统计量（F Stat）就是这个比例（用自由度修正过）。从上面结果看我们这个模型很显著（p值不超过万分之一），所以可以否定
第三类检验：
Type III Tests
Sum of Squares
Mean Square
这个表格给出了对各斜率项是否为零（ ）的检验结果。检验利用的是所谓第三类平方和（Type III SS），又叫偏平方和，它代表在只缺少了本变量的模型中加入本变量导致的模型平方和的增加量。比如，HEIGHT的第三类平方和即现在的模型平方和减去删除变量HEIGHT
的模型的模型平方和得到的差。第三类平方和与模型中自变量的次序无关，一般也不构成模型平方和的平方和分解。表中用F统计量对假设进行了检验，分子是第三类平方和的均方，分母为误差的均方。实际上，当分子自由度为1时，F统计量即通常的t检验统计量的平方。从表中可见，身高的作用是显著的，而年龄的作用则不显著，有可能去掉年龄后的模型更好一些。
参数估计及相关统计量：
Parameter Estimates
Parameter Estimates
Var Inflation
对截距项系数和各斜率项系数，给出了自由度（DF），估计值（Estimate），估计的标准误差（Std Error），检验系数为零的t统计量，t统计量的p值，检验共线性的容许度（Tolerance
）和方差膨胀因子（Var Inflation）。其中自变量 的容许度定义为1减去
对其它自变量的复相关系数平方，因此容许度越小（接近0），说明 对其它自变量的复相关系数平方大，即 可以很好地被其它自变量的线性组合近似，这样 在模型中的作用不大。记 ，则 ， 叫做方差膨胀因子，它代表 的系数估计的方差的比例系数，显然其值越大说明估计越不准确，也说明 在模型中的作用不大。方差膨胀因子与容许度互为倒数。
下一个结果为残差对预测值的散点图，用它可以检验残差中有无异常情况，比如非线性关系、异方差、模型辨识错误、异常值、序列相关等等。此例中各散点较随机地散布在0线的上下，没有明显的模式，可认为结果是合适的（多余的不显著的变量AGE不反映在残差图中）。
用Tables菜单可以加入一些其它的统计量。用Graphs菜单可以加入残差的正态概率图（Residual Normal QQ）和偏杠杆图（Partial
Leverage）。
在Vars菜单中可以指定一些变量，这些变量可以加入到数据窗口中。数据窗口的内容保存在内存中，不自动改写磁盘中的数据集，所以要保存数据窗口的修改结果的话需要用“File
| Save | Data”命令指定一个用来保存的数据集名。为了了解加入的变量的具体意义，选数据窗口菜单中的“Data Options”，选中“Show
Variable Labels”选项。各变量中，Hat Diag为帽子矩阵的对角线元素（帽子矩阵 恰好是 的）,即杠杆率，反映了每个观测的影响大小。Predicted为拟合值（预报值），Linear Predictor为使用线性模型拟合的结果，在线性回归时与Predicted相同。Residual为残差。Residual
Normal Quantile是残差由小到大排序后对应的标准正态的分位数，第 个残差的正态分位数用 计算，其中 为标准正态分布函数。Standardized
Residual （标准化误差）为残差除以其标准误差。Studentized Residual（学生化残差）为与标准化残差类似，但计算第 个学生化残差时预测值和方差估计都是在删除第 个观测后得到的。当学生化残差的值超过2时这个观测有可能是强影响点或异常点。
关于其它的一些诊断统计量请参考帮助菜单的“Extended Help | SAS System Help: Main menu | Help
for SAS Products | SAS/INSIGHT | Techniques | Multiple Regression” ，或《SAS系统：SAS/STAT软件使用手册》第一章和第九章。
在SAS/INSIGHT中，为了保存结果表格，在进行分析之前选中菜单“File | Save | Initial Tables”，这是一个状态开关，选中时输出表格画在分析窗口内的同时显示在输出（Output
）窗口。如果要保存某一个表格，也可以选定此表格（单击表格外框线），然后用菜单“File | Save | Tables”。为了保存分析窗口的图形，先选定此图形，然后选“File
| Save | Graphics File”，输入一个文件名，选择一种文件类型如BMP即可。为了打印某一表格或图形，先选定它，然后用菜单“File
| Print”。选中“File | Save | Statments”可以开始保存SAS/INSIGHT 语句。
经典线性回归理论的估计与假设检验要求自变量 为常数（非随机），随机误差项满足 。广义线性模型放宽了这些假设，其模型为
其中因变量 （ 向量）的元素为服从指数族分布（如正态、逆高斯、伽马、泊松、二项分布）的随机变量， （
向量）的元素为与 分布类型相同的随机误差项，元素之间相互独立，单调函数 叫做联系函数，它把因变量的均值 与自变量 （ 阵）的线性组合联系起来。 （ 向量）为回归系数。模型中每个自变量对应于设计阵 中的一列或几列， 的第一列一般元素全为1，对应于截距项。 （ 向量）是表示偏移量的变量。
注：随机变量Y称为服从指数族分布，如果其分布密度（概率函数）有如下形式：
其中 为自然参数或称经典参数， 为分散度参数（与尺度参数相关）， a, b,
c为确定性函数。这样的自变量Y的均值和方差与参数的关系如下：
为了使用SAS/INSIGHT拟合广义线性模型，在选“Analyze | Fit (Y X)”之后，选定因变量和自变量，然后按“Method”按钮，出现选择模型的对话框，在这里可以选因变量的分布类型（Response
Dist.），选联系函数，选估计尺度参数的方法。
各联系函数定义如下：
Probit ，其中 为标准正态分布函数
Comp. Log-Log
Power ， 在对话框的Power输入框指定。
对指数族中每一个因变量分布有一个特定的联系函数，使得 ，即用分布的期望值表示经典参数，这样的联系函数叫经典（canonical）联系函数。正态分布的经典联系函数为恒等变换，逆高斯分布为-2次方变换，伽玛分布为-1次方变换，泊松分布为对数变换，二项分布为逻辑变换（Logit
）。注意Logit、probit、复合重对数变换都只适用于二项分布。
例如，SASUSER.INGOTS中存放了一个铸造厂的数据，它记录了各批铸件在一定的加热、浸泡时间条件下出现的不能开始轧制的铸件数目。HEAT为加热时间，SOAK为浸泡时间，N为每批铸件的件数，R为加热浸泡后N件铸件中还不能开始轧制的铸件数。R应该服从二项分布，其分布参数（比例）可能受加热、浸泡时间的影响。因此，我们拟合以R为因变量，以HEAT和SOAK
为自变量的广义线性模型，因变量分布为二项分布，使用经典联系函数（Logit函数）。模型为
为了拟合这样的模型，选“Analyze | Fit(Y X)”，选R为Y变量，选HEAT和SOAK为自变量，按“Method”钮，选因变量分布为二项分布（Binomial），选变量N然后按“Binomial”钮，两次OK后即可以得到模型拟合窗口。可以看到，这个模型是显著的，但变量SOAK没有显著影响。去掉变量SOAK重新拟合模型。可以看出，HEAT的系数为0.0807是正数，说明加热时间越长不能轧制的件数越多。考察拟合结果窗口下方的残差对预报值图可以发现在右下方有三个异常点，用刷亮方法选定它们，可以看到，这三个观测都是总共只有一个铸件的，所以对一般结果意义不大。选“Edit
| Observations | Exclude in Calculation”可以把这几个点排除在外，发现结果基本不变。
SAS/STAT中提供了几个回归分析过程，包括REG（回归）、RSREG（二次响应面回归）、ORTHOREG（病态数据回归）、NLIN（非线性回归）、TRANSREG（变换回归）、CALIS（线性结构方程和路径分析）、GLM（一般线性模型）、GENMOD（广义线性模型），等等。我们这里只介绍REG过程，其它过程的使用请参考《SAS系统――SAS/STAT软件使用手册》。
REG过程的基本用法为：
PROC REG DATA=输入数据集 选项;
可参与建模的变量列表;
MODEL 因变量＝自变量表 / 选项;
REG过程是交互式过程，在使用了RUN语句提交了若干个过程步语句后可以继续写其它的REG 过程步语句，提交运行，直到提交QUIT语句或开始其它过程步或数据步才终止。
例如，我们对SASUSER.CLASS中的WEIGHT用HEIGHT和AGE建模，可以用如下的简单REG 过程调用：
proc reg data=sasuser.
model weight=
就可以在输出窗口产生如下结果，注意程序窗口的标题行显示“PROC REG Running”表示REG 过程还在运行，并没有终止。
Model: MODEL1
Dependent Variable: WEIGHT
Weight in pounds
Analysis of Variance
Parameter Estimates
Parameter=0
Prob & |T|
-141.223763
Height in inches
Age in years
这些结果与SAS/INSIGHT得到的结果是一致的。同样我们发现变量AGE的作用不显著，所以我们只要再提交如下语句：
model weight=
就可以得到第二个模型结果：
Model: MODEL2
Dependent Variable: WEIGHT
Weight in pounds
……………
事实上，REG提供了自动选择最优自变量子集的选项。在MODEL语句中加上“SELECTION= 选择方法”的选项就可以自动挑选自变量，选择方法有NONE（全用，这是缺省）、FORWARD
（逐步引入法）、BACKWARD（逐步剔除法）、STEPWISE（逐步筛选法）、MAXR（最大 增量法）、MINR（最小 增量法）、RSQUARE（ 选择法）、ADJRSQ（修正 选择法）、CP（Mallows的 统计量法）。比如，我们用如下程序：
model weight=height age / selection=
可得到如下结果：
Stepwise Procedure for Dependent Variable WEIGHT
Variable HEIGHT Entered
R-square = 0.
Sum of Squares
Mean Square
Regression
Sum of Squares
Bounds on condition number:
------------------------------------------------------------------------------
All variables left in the model are significant at the 0.1500 level.
No other variable met the 0.1500 significance level for entry into the model.
Summary of Stepwise Procedure for Dependent Variable WEIGHT
Entered Removed
Height in inches
可见只有变量HEIGHT进入了模型，而其它变量（AGE）则不能进入模型。
REG过程给出的缺省结果比较少。如果要输出高分辨率诊断图形的话需要在PROC REG 过程语句中加上GRAPHICS选项，用PRINT语句和PLOT语句显示额外的结果。为了显示模型的预测值（拟合值）和95％预测界限，使用语句
得到如下的结果：
Sum of Residuals
Sum of Squared Residuals
Predicted Resid SS (Press)
各列分别为观测序号（Obs），因变量的值（Dep Var），预测值（Predict Value），预测值的标准误差（Std Err Predict），95％预测区间下限（Lower
95% Predict），95％预测区间上限（Upper 95% Predict），残差（Residual，为因变量值减预测值）。在表后又给出了残差的总和（Sum
of Residuals），残差平方和（Sum of Squared Residuals），预测残差的平方和（Predicted Resid
SS (Press)）。所谓预测残差，是在计算第 i号观测的残差时从实际值中减去的预报值是用扣除第i号观测后的样本得到的模型产生的预报值，而不是我们一般所用的预测值（实际是拟合值）。第i号样本的预测残差还可以用公式
来计算，其中 为帽子矩阵 的第 i个主对角线元素。
用print cli列出的是实际值的预测界限，还可以列出模型均值的预测界限，使用
语句。在PRINT语句中可以指定的有ACOV, ALL, CLI, CLM, COLLIN, COLLINOINT, COOKD, CORRB,
COVB, DW, I, INFLUENCE, P, PARTIAL, PCORR1, PCORR2, R, SCORR1, SCORR2,
SEQB, SPEC, SS1, SS2, STB, TOL, VIF, XPX，等等。
对于自变量是一元的情况，可以在自变量和因变量的散点图上附加回归直线和均值置信界限。比如，
plot weight * height / conf95;
可以产生图 4，在图的上方列出了模型方程，右方还给出了观测个数、 、修正 、均方误差开根。在PLOT语句中可以使用PREDICTED. 、RESIDUAL.等特殊名字表示预测值、残差等计算出的变量，比如，在自变量为多元时无法作回归直线，常用的诊断图表为残差对预测值图，就可以用
plot residual. * predicted.;
绘制。为了绘制学生化残差的图形，可以用
plot rstudent. * obs.;
回归分析的其它用法及进一步的诊断方法请参考有关统计书籍和SAS使用手册。
统计学中用方差分析来研究分类变量（所谓“因素”）对数值型变量（所谓“指标” ）的影响。主要目的是研究某些因素对于指标有无显著的影响。对有显著影响的因素，一般希望找出最好水平。
单因素方差分析是4.1.2问题的一个自然延续。在4.1.2中，我们有一个分类变量把观测分为两组，我们要研究这两组的均值有没有显著差异。如果这个分类变量的取值不只两个，则这时4.1.2的检验方法不再适用，但我们同样要解决各组均值是否有显著差异的问题。如果各组之间有显著差异，说明这个因素（分类变量）对指标是有显著影响的，因素的不同取值（叫做水平）会影响到指标的取值。
例如，数据集SASUSER.VENEER中为比较若干种牌子的胶合板的耐磨情况得到的数据，变量BRAND 为试样的牌子，变量WEAR为试样的磨损量。共有五种牌子的胶合板，每种试验了4个试样。我们希望知道这五种牌子胶合板的磨损量有无显著差别，如果无显著差别我们在选购时就不必考虑哪一个更耐磨而只需考虑价格等因素，但如果结果有显著差异则应考虑使用耐磨性好的牌子。这里，因素是胶合板的牌子，指标为磨损量，当各种牌子胶合板磨损量有显著差异时，说明因素的取值对指标有显著的影响。所以，方差分析的结论是因素对指标有无显著影响。注意，经典的方差分析只判断因素的各水平有无显著差异，而不管两个因素之间是否有差异，比如说我们的五个牌子即使有四个牌子没有显著差异，只有一个牌子的胶合板比这四个牌子的都好，结论也是说因素是显著的，或因素的各水平间有显著差异。
方差分析把指标的方差分解为由因素的不同取值能够解释的部分，和剩余的不能解释的部分，然后比较两部分，当能用因素解释的部分明显大于剩余的部分时认为因素是显著的。方差分析假定观测是彼此独立的，观测为正态分布的样本，由因素各水平分成的各组的方差相等。在这些假定满足时，就可以用ANOVA过程来进行方差分析。其一般写法为
PROC ANOVA DATA＝数据集;
指标＝因素;
比如，为了分析SASUSER.VENEER中各种牌子的胶合板的耐磨性有无显著差别，首先我们假定假设检验使用的检验水平为0.05，可以使用如下程序进行方差分析：
proc anova data=sasuser.
model wear=
结果如下：
Analysis of Variance Procedure
Class Level Information
ACME AJAX CHAMP TUFFY XTRA
Number of observations in data set = 20
Dependent Variable: WEAR
Amount of material worn away
Sum of Squares
Mean Square
Corrected Total
Mean Square
结果可以分为四个部分，第一部分是因素水平的信息，我们看到因素只有一个BRAND，它有5个水平，分别是ACME、AJAX、CHAMP、TUFFY、XTRA。共有20个观测。第二部分就是经典的方差分析表，表前面指明了因变量（指标）为WEAR，第一列“来源”说明方差的来源，是模型的（可以用方差分析模型解释的），误差的（不能用模型解释的），还是总和。第三列为平方和，其大小代表了各方差来源作用的大小。第二列为自由度。第四列为均方，即平方和除以自由度。第五列F值是F统计量的值，其计算公式为模型均方除以误差均方，用来检验模型的显著性，如果不显著说明模型对指标的变化没有解释能力。第六列是F统计量的p值。由于这里p值小于0.05（我们的检验水平），所以模型是显著的，因素对指标有显著影响。结果的第三部分是一些与模型有关的简单统计量，第一个是复相关系数平方，与回归模型一样仍代表总变差中能被模型解释的比例，第二个是变异系数，第三个是根均方误差，第四个是指标的均值。结果的第四部分是方差分析表的细化，给出了各因素的平方和和F统计量，因为是单因素所以这一行与上面的“模型”一行相同。
当方差分析的正态分布假定或方差相等假定不能满足时，对单因素问题，可以使用非参数方差分析的Kruskal-Wallis检验方法。这种检验不要求观测来自正态分布总体，不要求各组的方差相等，甚至指标可以是有序变量（变量取值只有大小之分而没有差距的概念，比如磨损量可以分为大、中、小三档，得病的程度可以分为重、轻、无，等等）。
NPAR1WAY过程的调用与ANOVA过程不同，因为它是单因素方差分析过程，所以只要用CLASS 语句给出分类变量（因素），用VAR语句给出指标就可以了，一般格式为：
PROC NPAR1WAY DATA=数据集 WILCOXON;
注意这样的语句格式与1.4.2中两独立样本比较的做法完全相同。当“因素”有两个水平时，执行Wilcoxon秩和检验，多个水平时执行Kruskal-Wallis检验。
比如，为了分析上面的胶合板例子中各牌子的耐磨性有无显著差异，取定0.10的检验水平，可以用如下的NPAR1WAY过程：
proc npar1way data=sasuser.
得到如下结果：
N P A R 1 W A Y
P R O C E D U R E
Wilcoxon Scores (Rank Sums) for Variable WEAR
Classified by Variable BRAND
10.4830691
10.0000000
10.4830691
11.0000000
10.4830691
10.4830691
17.2500000
10.4830691
11.2500000
Average Scores Were Used for Ties
Kruskal-Wallis Test (Chi-Square Approximation)
Prob & CHISQ = 0.0175
结果分为两个部分，第一部分是各组的秩和的情况，包括观测个数（N）、秩和（Sum of Scores）、在各组无显著差异的零假设下的期望秩和（Expected
Under H0）、在零假设下的标准差（Std Dev Under H0）、平均秩和（Mean Score，为秩和除以组内观测数）。下面的“Average
Scores Were Used for Ties”是说当名次相同时（如两个第2）用名次的平均值（(2+3)/2=2.5）。第二部分为Kruskal-Wallis检验的结果，包括近似的
统计量，自由度，检验的p值（Prob & CHISQ ）。现在p值0.0175小于预定的水平0.10所以结论是各种牌子的胶合板的耐磨性能有显著差异。注意，Kruskal-Wallis检验是非参数检验，在同等条件下非参数检验一般比参数检验的功效低，所以这里的p值0.0175比用ANOVA过程得到的p值0.0017要大。
方差分析只检验各组是否没有任何两两之间的差异，但不检验到底是哪两组之间有显著差异。在三个或多个组之间进行两个或多个比较的检验叫做多重比较。多重比较在统计学中没有一个公认的解决方法，而是提供了若干种检验方法。因为多重比较要进行不只一次的比较，所以在多重比较的检验水平有两种：总错误率（experimentwise
error rate）和单次比较错误率。总错误率是指所有比较（比如，五个组两两之间比较有10次）的总第一类错误概率，单次比较错误率是指每一次比较的第一类错误概率。显然，总错误率要比单次比较错误率高。
在ANOVA过程中使用MEANS语句可以进行多重比较。格式如下：
因素 / 选项;
如果不使用选项，则只对因素的各水平计算指标的平均值和标准差，比如：
proc anova data=sasuser.
model wear=
则在通常的方差分析结果基础上增加如下结果：
-------------WEAR------------
要进行两两比较，有多种方法，可以在MEANS语句的选项中指定检验方法。
一、用重复t检验控制单次比较错误率
重复t检验的想法很简单：在适当的检验水平下对两组之间进行两样本t检验并对所有组两两之间检验。控制的是每次比较的第一类错误概率。缺省使用0.05水平。注意这样检验的总错误率将大大高于每次比较的错误率。比如，在上面程序后加入（ANOVA是交互式过程）
means brand /
可得如下结果：
T tests (LSD) for variable: WEAR
NOTE: This test controls the type I comparisonwise error rate not the
experimentwise error rate.
Alpha= 0.05
MSE= 0.020833
Critical Value of T= 2.13
Least Significant Difference= 0.2175
Means with the same letter are not significantly different.
T Grouping
结果先说明了检验的指标是变量WEAR，然后说明了这种检验控制单次比较的第一类错误概率而不是总的第一类错误概率。下面给出了检验的一些指标，比如水平（Alpha）为0.05（控制单次比较的第一类错误概率），自由度（df）为15，误差的均方（MSE，是方差分析表中误差的均方）为0.020833，两样本t检验的t统计量的临界值（Critical
Value of T）为2.13，如果两样本t检验的t统计量值绝对值超过临界值则认为两组有显著差异，或者等价地，如果两组的均值之差绝对值大于最小显著差别（Least
Significant Difference）0.2175也是有显著差异。下面列出了检验的结果，把因素各水平的指标平均值由大到小排列，然后把两两比较的结果用第一列的字母来表示，字母相同的水平没有显著差异，字母不同的水平有显著差异。所以我们看到，重复t检验的结果把五种牌子分成了A、B、C三个组，TUFFY单独是一组，它的磨损量最大；XTRA、CHAMP、ACME是一组，这三种两两之间没有显著差异；AJAX单独是一组，其磨损量最小。
二、用Bonferroni t检验控制总错误率
Bonferroni t检验通过把每次比较的错误率取得很小来控制总误差率。比如，共有10次比较时，把每次比较的错误率控制在0.005就可以保证总错误率不超过0.05，但是，这样得到的实际总第一类错误率可能要比预定的水平小得多。在MEANS语句中使用BON语句可以执行Bonferroni
t检验，缺省总错误率控制水平为0.05。对上面的胶合板数据增加如下语句：
means brand /
结果如下：
Bonferroni (Dunn) T tests for variable: WEAR
NOTE: This test controls the type I experimentwise error rate, but
generally has a higher type II error rate than REGWQ.
Alpha= 0.05
MSE= 0.020833
Critical Value of T= 3.29
Minimum Significant Difference= 0.3354
Means with the same letter are not significantly different.
Bon Grouping
结果先说明了检验类型和指标（变量WEAR），然后说明了检验控制总第一类错误率，但一般比REGWQ方法的第二类错误概率高（检验功效较低）。下面给出了几个检验用的值。最后给出了Bonferroni
t检验的结果，有相同分组字母的因素水平间无显著差异，否则有显著差异。我们看到，TUFFY与XTRA、CHAMP、ACME没有显著差异，与AJAX有显著差异；XTRA、CHAMP、ACME
两两之间没有显著差异，而且与其它两个也都没有显著差异；AJAX与TUFFY有显著差异，与其它三个没有显著差异。其分组是有交叉的。
三、用REGWQ检验控制总错误率
用Bonferroni t检验控制总错误率过于保守，功效较低，不易发现实际存在的显著差异。REGWQ方法可以控制总错误率并且一般比Bonferroni
t检验要好。这种方法执行多阶段的检验，它对因素水平的各种子集进行检验。在MEANS语句中用REGWQ选项可以进行REGWQ检验。例如，在前面的例子后再运行
means brand/
Ryan-Einot-Gabriel-Welsch Multiple Range Test for variable: WEAR
NOTE: This test controls the type I experimentwise error rate.
Alpha= 0.05
MSE= 0.020833
Number of Means
Critical Range
0...2941581
Means with the same letter are not significantly different.
REGWQ Grouping
可见它比Bonferroni方法发现了较多的显著差异，除了TUFFY和AJAX仍有显著差异以外，还发现XTRA、CHAMP、ACME也都与AJAX有显著差异。
MEANS语句的选项可以同时使用。在MEANS语句中可以用ALPHA=水平值来指定检验的水平。ANOVA过程中还提供了其它的多重比较方法，请自己参考有关资料。
SAS提供了若干个方差分析过程，可以考虑多个因素、有交互作用、有嵌套等情况的方差分析。用GLM过程还可以用一般线性模型来处理方差分析问题。在这里我们只介绍如何用ANOVA
过程进行均衡设计的多因素方差分析。
例如，为了提高一种橡胶的定强，考虑三种不同的促进剂（因素A）、四种不同分量的氧化锌（因素B）对定强的影响，对配方的每种组合重复试验两次，总共试验了24次，得到如下的结果：
表格 1 橡胶配方试验数据
我们首先把数据输入为SAS数据集。输入的办法可以是直接输入各个观测，例如：
input A B STREN;
……………
也可以使用如下的直接循环控制的INPUT读取：
input stren @@;
其中INPUT语句尾部的两个@符号表示多次INPUT语句可以从同一行去读取（否则每次INPUT 语句运行时自动从下一行开始读）。
为了研究两个因素的主效应和交互作用，使用如下ANOVA过程：
proc anova data=
stren = a b a*b;
MODEL语句中中A表示因素A的主效应，B表示因素B的主效应，A*B表示A和B的交互作用。结果如下：
Analysis of Variance Procedure
Class Level Information
Number of observations in data set = 24
Dependent Variable: STREN
Corrected Total
STREN Mean
Mean Square
结果首先给出了因素（Class）的变量名和各水平值，观测数。然后是总的方差分析表，指明指标为变量STREN，给出了模型、误差、总平方和，F统计量值和p值。可见模型是显著的。为了分析各作用的显著性，看后面的详细的方差分析表，它给出了模型中各作用（A、B、A*B
）的平方和和检验的F统计量值及p值。可以看出，两个因素的主效应都是显著的，交互作用效应不显著。所以，我们可以重新运行ANOVA过程，不指定交互作用效应：
proc anova data=
stren = a b /
这时模型的F统计量变为30.53，因素A主效应的F统计量变为22.89，因素B主效应的F统计量变为35.63，都增大了。两个因素的主效应仍是高度显著的，说明它们对定强都有显著影响。为了找到最好的配方，在前面的ANOVA过程后使用
可以计算出每种水平下的指标平均值，因素A（促进剂）在第三水平使指标（定强）最大，因素B（氧化锌）在第四水平使指标最大，所以最好的配方是：第三种促进剂，第四种氧化锌分量。
ANOVA也可以用来分析正交设计的结果。例如，为了提高某种试剂产品的收率（指标），考虑如下几个因素对其影响：
表格 2 试剂产品影响因素
A：反应温度
B：反应时间
C：硫酸浓度
D：硫酸产地
E：操作方式
2 (不搅拌)
把这五个因素放在 表的五列上，得到如下的试验方案及结果（见下面的数据步）。用ANOVA过程可以分析：
input temp time
proc anova data=
model prod = temp--
means temp--mix /
用0.05水平，得到的模型是显著的（模型p值为0.0250），各因素的检验结果如下：
Mean Square
可见硫酸浓度、产地、操作方式是显著的，必须采用它们的最好水平，温度、时间不显著，在同等条件下可以优先采用它们的最好水平。从MEANS语句的结果可以知道，硫酸浓度的最好水平是水平2（27%），硫酸产地的最好水平是水平2（上海），操作方式的最好水平是水平2
（不搅拌），反应温度的最好水平是水平1（50℃），反应时间的最好水平是水平1（1小时）。从以上分析可以得到好的生产方案。
上面所讲到的统计分析主要针对数值型（区间）变量进行。在实际工作中，离散取值的名义变量（如性别、职业、民族）和有序变量（如调查意见的完全同意、同意、中立、反对、强烈反对）也是十分常见的，对这类离散变量（又称属性变量）的分析也是统计学的重要的研究内容。这一节我们讲述检验两个离散取值的变量的独立性的列联表检验方法，并介绍有序变量的关联性量度的算法。
离散变量的取值可以把样本进行分类。比如，我们的样本是一个班的学生情况，可以根据学生的性别把观测分为男生和女生两个组。我们也可以根据学生的来源把观测分为本地学生和外地学生两个组。如果联合使用这两个变量对观测分类，就可以把观测分为四个组，我们可以统计每个组学生的人数，并把结果画成一个表格：
学生性别、来源分布表
这样的表格就叫做列联表。它给出了按照两个变量总和分类得到的每一个小类的观测个数。
为了得到这样的表格，需要把数据输入为数据集。有时我们得到的数据是每一个观测的变量取值，比如，我们有每一个学生的性别（SEX）情况和来源（FROM）情况，可以输入这些原始数据，如：
input sno sex $ from $;
label sex='性别' from='来源';
…………/* 所有学生的记录 */
然后用如下的FREQ过程可以画出列联表：
proc freq data=
tables from *
结果见表格4。
表格数据输入的另一种情况是，我们得到的数据就已经是上面表格 3那样的调查结果而不是具体的样本情况，可以直接把表格输入一个数据集，但数据集中要有一个代表观测数的变量，例如：
input from $ sex $
label sex='性别' from='来源';
这样的数据要画列联表，需要在FREQ过程中使用WEIGHT语句指定表示重复数的变量（NUMCELL ）：
proc freq data=
tables from *
结果和上面得到的结果相同。
在输出结果中，我们看到TABLES语句中的前一个变量被用来区分行，后一个变量被用来区分列。每个格子中有四个数：Frequency（频数，本格子的观测数），Percent（百分比），Row
Pct（行百分比，表示本类在本行中占的百分比，比如本地男生有4个人，本行有10 个人，占本行的40.00%），Col Pct（列百分比）。在表的右侧有行总计，比如本地学生有10
个人，占总学生数（31人）的32.26%。在表的下侧有列总计，比如男生有18个人，占学生总数的58.06%。表格右下方是总数（31）和总百分比（100）。
为了作列联表，调用FREQ过程，使用TABLES语句指定行变量和列变量，两者用星号分开，如果数据本身是表格数据还需要用WEIGHT语句指定存放表格单元观测数的变量。
可以作出简化的表格，在TABLES语句中加上NOFREQ、NOPCT、NOROW、NOCOL等选项就可以抑制相应的统计量的输出。例如，用如下程序：
proc freq data=
tables from * sex /
就可以产生只有单元数的表格。
对于数值型变量，我们考虑其相关关系的通常的办法是计算相关系数和进行回归分析。如果我们要研究离散取值的名义变量和有序变量有无相关，最常用的检验办法是列联表独立性检验。列联表检验的零假设是两变量
和 相互独立，计算一个 统计量，与列联表中频数取值和零假设下期望取值之差有关，当 很大时否定零假设。
例如，为了探讨吸烟与慢性支气管炎有无关系，调查了339人，情况如下：
吸烟与慢性支气管炎调查表
患慢性支气管炎
未患慢性支气管炎
设想有两个随机变量X，Y：X取1表示吸烟，取2表示不吸烟，Y取1表示患慢性支气管炎，取2表示未患。零假设为：
: X与Y相互独立
要检验此零假设，先取定检验水平0.05，用PROC FREQ过程，在TABLES语句中加上CHISQ选项即可。下面的例子中还加入了EXPECTED选项以显示零假设下的期望频数值：
input smoke $ bron $
label smoke='吸烟'
bron='慢性支气管炎';
proc freq data=
tables smoke*bron / nopct norow n
结果如下：
TABLE OF SMOKE BY BRON
SMOKE(吸烟)
BRON(慢性支气管炎)
Frequency|
Expected |患病
---------+--------+--------+
| 22.136 | 111.86 |
---------+--------+--------+
| 33.864 | 171.14 |
---------+--------+--------+
STATISTICS FOR TABLE OF SMOKE BY BRON
------------------------------------------------------
Chi-Square
Likelihood Ratio Chi-Square
Continuity Adj. Chi-Square
Mantel-Haenszel Chi-Square
Fisher's Exact Test (Left)
Phi Coefficient
Contingency Coefficient
Cramer's V
Sample Size = 339
列联表中列出了表格单元频数和在零假设下的期望频数，可以看出，吸烟人中患病的数目比期望数目大。检验的结果只要看后面的统计量部分的Chi-Square一行，其值为7.469，p值为0.006，所以应否定零假设，吸烟与患慢性支气管炎是不独立的。
使用 检验要求每个单元格至少频数不少于5。在条件不满足的时候还可以使用Fisher精确检验。对于两行两列的表格FREQ过程自动给出Fisher
精确检验的结果，其双侧检验p值为0.00686，应拒绝零假设。
对于区间变量，我们可以计算两两的相关系数。属性变量因为没有数值概念所以不能计算相关系数，但对于两个有序变量我们可以计算类似于相关系数的关联性量度。其中一种关联性量度叫做Kendal
Tau-b统计量，取值在-1到1之间，值接近于1表示正关联，接近于-1 表示负关联，接近于0表示没有相关关系。
下面用例子说明如何在FREQ过程中计算Kendal Tau-b统计量。本例取自《SAS系统与基础统计分析》一书。假设我们要研究奶牛种群大小与其患某种细菌性疾病的关系。牛的患病程度（DISEASE）分为没有（0）、低（1）、高（2），牛群大小（HERDSIZE）分为小（1）、中（2）、大（3）。数据如下数据步所示：
input herdsi
label herdsize='牛群大小'
disease='患病程度';
用FREQ过程在TABLES语句中加上MEASURES选项就可以计算Kendall Tau-b统计量：
proc freq data=
tables herdsize*disease / measures expected n
title '奶牛疾病数据分析';
结果如下：
HERDSIZE(牛群大小)
DISEASE(患病程度)
Frequency|
Expected |
---------+--------+--------+--------+
| 2.9231 | 7.4615 | 12.615 |
---------+--------+--------+--------+
| 5.2107 | 13.301 | 22.488 |
---------+--------+--------+--------+
| 29.866 | 76.237 |
---------+--------+--------+--------+
STATISTICS FOR TABLE OF HERDSIZE BY DISEASE
------------------------------------------------------
Kendall's Tau-b
Stuart's Tau-c
Somers' D C|R
Somers' D R|C
Pearson Correlation
Spearman Correlation
Lambda Asymmetric C|R
Lambda Asymmetric R|C
Lambda Symmetric
Uncertainty Coefficient C|R
Uncertainty Coefficient R|C
Uncertainty Coefficient Symmetric
Sample Size = 299
算出的Kendall Tau-b统计量值为0.217，渐近标准误差（ASE）为0.061，用统计量值加减两倍标准误差作为Kendall Tau-b的95%置信区间，可算得(0.095，0.339)在零点左边，所以可认为奶牛患病程度与种群大小有正的关联。事实上，我们从列联表中实际频数与期望频数的对比也可以看出，小的种群患病比期望值轻，大的种群患病比期望值重，即患病程度与种群大小有正的关联。