统计学的统一

Posted 2022-01-04 sssr

很久之间就看过https://bbs.pinggu.org/thread-2736449-1-1.html 这篇帖子，感觉很有意思，最近正好研究广义估计方程，就把这套理论用sas实现以下。

受到该启发后，认真再复习GLM的相关资料，得到更加重要的总结如下（来自高惠璇SAS/STAT软件使用手册，实际是SAS8.2的User‘s guide的中文版，但是目前SAS 9.2，9.3的User‘s guide关于GLM模型的介绍中已经删去了这么经典的总结，实在可惜，倒让人看不到GLM的真正长处了）：
如果X1-X3，Y1-Y2为连续性变量，Y3为分类变量，a-c为分类变量，time为时间变量，目前我们熟悉的模型可以简单概括如下：
(1) y1= x1    简单回归
(2) y1= x1 x2 x3  多重回归（multiple regression）
(3) y1 y2=x1 x2   多元回归（multivariate regression）
(4) y1= a   单因素方差分析
(5) y1= a b   （析因设计的）主效应分析
(6) y1= a b  a*b （析因设计的）主效应加交互项分析
(7) y1= a x1   协方差分析
(8) y3= a   单因素logistic回归
(9) y3= a b c x1 x2 x3   多因素logistic回归
(10) y3(time) =a   单因素cox回归
(11) y3(time) = a b c x1 x2 x3  多因素cox回归

1-7采用SAS的一般线性模型GLM都能实现，而1-9采用SAS的广义线性模型GENMOD都能实现。

转载请注明出处：https://www.cnblogs.com/SSSR/p/11006749.html 

以下为两组的t检验（没做正态性和方差齐性），为例
data a;
input group y;
cards;
1 11.1
1 8.17
1 12.73
1 15.83
1 15.6
1 17.2
1 5.45
1 11.3
2 23.73
2 18.86
2 26.65
2 16.72
2 17.33
2 18.08
2 16.55
2 17.87
;
run;
PROC TTEST DATA=a;
TITLE "T-test Example";
CLASS GROUP;
VAR y;
RUN;

 

方法	方差	自由度	t 值	Pr > |t|
汇总	等于	14	-3.78	0.0020
Satterthwaite	不等于	13.89	-3.78	0.0021

proc anova data=a;
class group;
model y=group;
means group;
means group/hovtest;
run;
quit;

源	自由度	Anova 平方和	均方	F 值	Pr > F
group	1	213.2330063	213.2330063	14.28	0.0020

F值=t值的平方

proc reg;/*调用回归模块*/
model y=group/cli; /*对y关于x1做回归，/cli表示求预测值与预测区间*/
run;
quit;

 

参数估计
变量	自由度	参数 估计	标准 误差	t 值	Pr > |t|
Intercept	1	4.87125	3.05506	1.59	0.1331
group	1	7.30125	1.93219	3.78	0.0020

proc genmod data=a;
class group(ref=first);
model y=group;
run;

 

最大似然参数估计的分析
参数		自由度	估计	标准 误差	Wald 95% 置信限		Wald 卡方	Pr > 卡方
Intercept		1	12.1725	1.2780	9.6676	14.6774	90.72	<.0001
group	2	1	7.3013	1.8074	3.7588	10.8437	16.32	<.0001
group	1	0	0.0000	0.0000	0.0000	0.0000	.	.
尺度		1	3.6148	0.6390	2.5563	5.1116

这里参数估计是一样的，但是用的是卡方做检验，不知道卡方值和T值 F值之间有没有什么关系？？？？？

proc gam data=a;
model y=param(group);
run;

gam是线性可加模型，用来判断自变量与因变量之间是否存在非线性关系的，这里使用param则表示不进行多项式判断。得到的结果reg一致。

回归模型分析 参数估计
参数	参数 估计	标准 误差	t 值	Pr > |t|
Intercept	4.87125	3.05506	1.59	0.1331
group	7.30125	1.93219	3.78	0.0020

 

 

再用一个例子，对比一下gam、genmod和reg之间的关系

 

title ‘Patterns of Diabetes‘;
data diabetes;
input Age BaseDeficit CPeptide @@;
logCP = log(CPeptide);
datalines;
5.2 -8.1 4.8 8.8 -16.1 4.1 10.5 -0.9 5.2
10.6 -7.8 5.5 10.4 -29.0 5.0 1.8 -19.2 3.4
12.7 -18.9 3.4 15.6 -10.6 4.9 5.8 -2.8 5.6
1.9 -25.0 3.7 2.2 -3.1 3.9 4.8 -7.8 4.5
7.9 -13.9 4.8 5.2 -4.5 4.9 0.9 -11.6 3.0
11.8 -2.1 4.6 7.9 -2.0 4.8 11.5 -9.0 5.5
10.6 -11.2 4.5 8.5 -0.2 5.3 11.1 -6.1 4.7
12.8 -1.0 6.6 11.3 -3.6 5.1 1.0 -8.2 3.9
14.5 -0.5 5.7 11.9 -2.0 5.1 8.1 -1.6 5.2
13.8 -11.9 3.7 15.5 -0.7 4.9 9.8 -1.2 4.8
11.0 -14.3 4.4 12.4 -0.8 5.2 11.1 -16.8 5.1
5.1 -5.1 4.6 4.8 -9.5 3.9 4.2 -17.0 5.1
6.9 -3.3 5.1 13.2 -0.7 6.0 9.9 -3.3 4.9
12.5 -13.6 4.1 13.2 -1.9 4.6 8.9 -10.0 4.9
10.8 -13.5 5.1
;

proc gam data=diabetes;
model logCP = spline(Age) spline(BaseDeficit);
run;

 

回归模型分析 参数估计
参数	参数 估计	标准 误差	t 值	Pr > |t|
Intercept	1.48141	0.05120	28.93	<.0001
Linear(Age)	0.01437	0.00437	3.28	0.0024
Linear(BaseDeficit)	0.00807	0.00247	3.27	0.0025

proc genmod data=diabetes;
model logcp=age basedeficit;
run;

 

最大似然参数估计的分析
参数	自由度	估计	标准 误差	Wald 95% 置信限		Wald 卡方	Pr > 卡方
Intercept	1	1.4829	0.0576	1.3700	1.5957	663.41	<.0001
Age	1	0.0150	0.0049	0.0054	0.0247	9.32	0.0023
BaseDeficit	1	0.0090	0.0028	0.0035	0.0144	10.42	0.0012
尺度	1	0.1247	0.0134	0.1009	0.1540

proc reg data=diabetes;
model logcp=age basedeficit;
run;

 

大家可以看到三者的结果是一致的。知识genmod的检验方法不同而已。

参数估计
变量	自由度	参数 估计	标准 误差	t 值	Pr > |t|
Intercept	1	1.48285	0.05969	24.84	<.0001
Age	1	0.01502	0.00510	2.95	0.0054
BaseDeficit	1	0.00896	0.00288	3.11	0.0034