《R语言》第9章：回归分析

返回课程主页

版权声明：本套课程材料开源，使用和分享必须遵守「创作共用许可协议 CC BY-NC-SA」（来源引用-非商业用途使用-以相同方式共享）。

---

Chap09：回归分析

往期要点回顾

• Chap08 # 多因素方差分析
• Chap08 # 简单效应检验与多重比较

本章要点目录

• 【知识点】回归分析（regression）统计术语
• 【实践1】截距模型
• 【实践2】一元回归模型
• 【实践3】多元回归模型
• 【实践4】交互作用模型（重点）
• 【知识点】为什么回归分析交互作用中的斜率并不是“主效应”？（重点）
• 【实践5】如何让回归分析交互作用中的斜率具有“主效应”含义？
• 【实践6】多元回归模型的结果表输出（重点）
• 【实践7】空气质量数据回归分析实例
• 【实践8】基于PROCESS()函数的调节效应与简单斜率分析（重点）
• 【知识点】回归模型的基本假定
• 【实践9】检查回归模型是否满足基本假定
• 【探索发现】各类回归模型使用的R函数

## 本章所需R包
library(bruceR)  # 回归模型结果的 Word 表格输出
library(sjPlot)  # 回归模型结果的 HTML 表格输出

回归分析的统计知识回顾

【知识点】回归分析（regression）统计术语

多元回归模型的一般形式：

\[ Y_i = b_0 + b_1 X_1 + b_2 X_2 + ... + b_p X_p + e_i \]

• 截距（intercept）
  • \(b_0\)：所有自变量都取0时，因变量的预测期望值
    • 如果没有任何自变量，截距就是因变量的平均值
• 斜率（slope）
  • \(b_1\), …, \(b_p\)：自变量每变化1，因变量变化多少个单位
    • 非标准化系数：\(b\)
    • 标准化系数：\(\beta\)
• 残差（residual）
  • \(e_i = Y_i - \hat{Y_i}\)：每个个体（\(i\)）因变量的真实观测值（\(Y_i\)）与模型预测值（\(\hat{Y_i}\)）之间的差异

回归模型建立与结果报告

## 数据准备
data = data.table(
  A = c(0, 0, 0,  0, 0, 0,  1, 1, 1,  1, 1, 1),
  B = c(0, 0, 0,  1, 1, 1,  0, 0, 0,  1, 1, 1),
  Y = c(0, 1, 2,  3, 4, 5,  0, 1, 2,  6, 7, 8)
)
data

        A     B     Y
    <num> <num> <num>
 1:     0     0     0
 2:     0     0     1
 3:     0     0     2
 4:     0     1     3
 5:     0     1     4
 6:     0     1     5
 7:     1     0     0
 8:     1     0     1
 9:     1     0     2
10:     1     1     6
11:     1     1     7
12:     1     1     8

## 方差分析
m = MANOVA(data, dv="Y", between=c("A", "B"))

====== ANOVA (Between-Subjects Design) ======

Descriptives:
────────────────────────
 "A" "B"  Mean    S.D. n
────────────────────────
  A0  B0 1.000 (1.000) 3
  A0  B1 4.000 (1.000) 3
  A1  B0 1.000 (1.000) 3
  A1  B1 7.000 (1.000) 3
────────────────────────
Total sample size: N = 12

ANOVA Table:
Dependent variable(s):      Y
Between-subjects factor(s): A, B
Within-subjects factor(s):  –
Covariate(s):               –
─────────────────────────────────────────────────────────────────────
           MS   MSE df1 df2      F     p     η²p [90% CI of η²p]  η²G
─────────────────────────────────────────────────────────────────────
A       6.750 1.000   1   8  6.750  .032 *     .458 [.034, .715] .458
B      60.750 1.000   1   8 60.750 <.001 ***   .884 [.680, .940] .884
A * B   6.750 1.000   1   8  6.750  .032 *     .458 [.034, .715] .458
─────────────────────────────────────────────────────────────────────
MSE = mean square error (the residual variance of the linear model)
η²p = partial eta-squared = SS / (SS + SSE) = F * df1 / (F * df1 + df2)
ω²p = partial omega-squared = (F - 1) * df1 / (F * df1 + df2 + 1)
η²G = generalized eta-squared (see Olejnik & Algina, 2003)
Cohen’s f² = η²p / (1 - η²p)

Levene’s Test for Homogeneity of Variance:
───────────────────────────────────
       Levene’s F df1 df2     p    
───────────────────────────────────
DV: Y       0.000   3   8 1.000    
───────────────────────────────────

emmip(m, ~ A, CIs=TRUE)

NOTE: Results may be misleading due to involvement in interactions

emmip(m, ~ B, CIs=TRUE)

NOTE: Results may be misleading due to involvement in interactions

emmip(m, ~ A | B, CIs=TRUE)

emmip(m, B ~ A, CIs=TRUE)

【实践1】截距模型

## lm(formula=..., data=...)

m1 = lm(Y ~ 1, data)  # 1表示模型截距
summary(m1)

Call:
lm(formula = Y ~ 1, data = data)

Residuals:
   Min     1Q Median     3Q    Max 
 -3.25  -2.25  -0.75   2.00   4.75 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)   
(Intercept)    3.250      0.789    4.12   0.0017 **
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.73 on 11 degrees of freedom

bruceR::model_summary(m1)

Model Summary

──────────────────────
             (1) Y    
──────────────────────
(Intercept)   3.250 **
             (0.789)  
──────────────────────
R^2           0.000   
Adj. R^2      0.000   
Num. obs.    12       
──────────────────────
Note. * p < .05, ** p < .01, *** p < .001.

Not enough model terms in the conditional part of the model to check for
  multicollinearity.

sjPlot::tab_model(m1)

	Y
Predictors	Estimates	CI	p
(Intercept)	3.25	1.51 – 4.99	0.002
Observations	12
R2 / R2 adjusted	0.000 / 0.000

sjPlot::tab_model(
  m1, digits=3, show.ci=FALSE,
  show.se=TRUE, collapse.se=TRUE
)

	Y
Predictors	Estimates	p
(Intercept)	3.250 (0.789)	0.002
Observations	12
R2 / R2 adjusted	0.000 / 0.000

## 截距含义：Y的均值
mean(data$Y)  # 截距等于均值

[1] 3.25

predict(m1)  # 预测期望值

   1    2    3    4    5    6    7    8    9   10   11   12 
3.25 3.25 3.25 3.25 3.25 3.25 3.25 3.25 3.25 3.25 3.25 3.25 

TTEST(data, y="Y", digits=3)  # 截距模型等价于单样本t检验

One-Sample t-test

Hypothesis: two-sided (μ ≠ 0)

Descriptives:
──────────────────────────
 Variable  N   Mean (S.D.)
──────────────────────────
        Y 12 3.250 (2.734)
──────────────────────────

Results of t-test:
────────────────────────────────────────────────────────────────────────
                t df     p      Difference [95% CI]   Cohen’s d [95% CI]
────────────────────────────────────────────────────────────────────────
Y: (Y - 0)  4.117 11  .002 **  3.250 [1.513, 4.987] 1.189 [0.553, 1.824]
────────────────────────────────────────────────────────────────────────

【实践2】一元回归模型

m2 = lm(Y ~ 1 + A, data)
m2 = lm(Y ~ A, data)  # 相当于：Y ~ 1 + A
# 当有自变量时，截距1默认存在，可省略不写
summary(m2)

Call:
lm(formula = Y ~ A, data = data)

Residuals:
   Min     1Q Median     3Q    Max 
 -4.00  -2.12   0.00   2.12   4.00 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)  
(Intercept)     2.50       1.12    2.23     0.05 *
A               1.50       1.59    0.95     0.37  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.75 on 10 degrees of freedom
Multiple R-squared:  0.0821,    Adjusted R-squared:  -0.00973 
F-statistic: 0.894 on 1 and 10 DF,  p-value: 0.367

tab_model(m2, digits=3, show.ci=FALSE,
          show.se=TRUE, collapse.se=TRUE)

	Y
Predictors	Estimates	p
(Intercept)	2.500 (1.122)	0.050
A	1.500 (1.586)	0.367
Observations	12
R2 / R2 adjusted	0.082 / -0.010

## 截距含义：A=0时，Y的预测期望值（等于均值）
mean(data[A==0]$Y)

[1] 2.5

predict(m2)  # 预测期望值

  1   2   3   4   5   6   7   8   9  10  11  12 
2.5 2.5 2.5 2.5 2.5 2.5 4.0 4.0 4.0 4.0 4.0 4.0 

predict(m2, data.table(A=0))  # A=0时，Y的预测期望值

  1 
2.5 

## 斜率含义：A从0变到1，Y的均值变化
TTEST(data, y="Y", x="A", digits=3)  # 一元回归等价于独立样本t检验

Independent-Samples t-test

Hypothesis: two-sided (μ2 - μ1 ≠ 0)

Descriptives:
──────────────────────────────────────
 Variable Factor Level N   Mean (S.D.)
──────────────────────────────────────
        Y      A     0 6 2.500 (1.871)
        Y      A     1 6 4.000 (3.406)
──────────────────────────────────────

Levene’s test for homogeneity of variance:
──────────────────────────────────────────
              Levene’s F df1 df2     p    
──────────────────────────────────────────
Y: A (1 - 0)        8.44   1  10  .016 *  
──────────────────────────────────────────
Note: H0 = equal variance (homoscedasticity).
If significant (violation of the assumption),
then you should better set `var.equal=FALSE`.

Results of t-test:
────────────────────────────────────────────────────────────────────────────
                  t df     p       Difference [95% CI]    Cohen’s d [95% CI]
────────────────────────────────────────────────────────────────────────────
Y: A (1 - 0)  0.946 10  .367     1.500 [-2.035, 5.035] 0.546 [-0.741, 1.832]
────────────────────────────────────────────────────────────────────────────

m3 = lm(Y ~ 1 + B, data)
m3 = lm(Y ~ B, data)  # 相当于：Y ~ 1 + B
# 当有自变量时，截距1默认存在，可省略不写
summary(m3)

Call:
lm(formula = Y ~ B, data = data)

Residuals:
   Min     1Q Median     3Q    Max 
  -2.5   -1.0    0.0    1.0    2.5 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)    1.000      0.599    1.67  0.12577    
B              4.500      0.847    5.32  0.00034 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.47 on 10 degrees of freedom
Multiple R-squared:  0.739, Adjusted R-squared:  0.712 
F-statistic: 28.3 on 1 and 10 DF,  p-value: 0.00034

tab_model(m3, digits=3, show.ci=FALSE,
          show.se=TRUE, collapse.se=TRUE)

	Y
Predictors	Estimates	p
(Intercept)	1.000 (0.599)	0.126
B	4.500 (0.847)	<0.001
Observations	12
R2 / R2 adjusted	0.739 / 0.712

## 截距含义：B=0时，Y的预测期望值（等于均值）
mean(data[B==0]$Y)

[1] 1

predict(m3)  # 预测期望值

  1   2   3   4   5   6   7   8   9  10  11  12 
1.0 1.0 1.0 5.5 5.5 5.5 1.0 1.0 1.0 5.5 5.5 5.5 

predict(m3, data.table(B=0))  # B=0时，Y的预测期望值

1 
1 

## 斜率含义：B从0变到1，Y的均值变化
TTEST(data, y="Y", x="B", digits=3)  # 一元回归等价于独立样本t检验

Independent-Samples t-test

Hypothesis: two-sided (μ2 - μ1 ≠ 0)

Descriptives:
──────────────────────────────────────
 Variable Factor Level N   Mean (S.D.)
──────────────────────────────────────
        Y      B     0 6 1.000 (0.894)
        Y      B     1 6 5.500 (1.871)
──────────────────────────────────────

Levene’s test for homogeneity of variance:
──────────────────────────────────────────
              Levene’s F df1 df2     p    
──────────────────────────────────────────
Y: B (1 - 0)        3.91   1  10  .076 .  
──────────────────────────────────────────
Note: H0 = equal variance (homoscedasticity).
If significant (violation of the assumption),
then you should better set `var.equal=FALSE`.

Results of t-test:
──────────────────────────────────────────────────────────────────────────
                  t df     p      Difference [95% CI]   Cohen’s d [95% CI]
──────────────────────────────────────────────────────────────────────────
Y: B (1 - 0)  5.316 10 <.001 *** 4.500 [2.614, 6.386] 3.069 [1.783, 4.355]
──────────────────────────────────────────────────────────────────────────

【实践3】多元回归模型

m4 = lm(Y ~ A + B, data)  # 相当于：Y ~ 1 + A + B
# 当有自变量时，截距1默认存在，可省略不写
summary(m4)

Call:
lm(formula = Y ~ A + B, data = data)

Residuals:
   Min     1Q Median     3Q    Max 
 -1.75  -0.75   0.00   0.75   1.75 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)    0.250      0.640    0.39  0.70520    
A              1.500      0.739    2.03  0.07300 .  
B              4.500      0.739    6.09  0.00018 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.28 on 9 degrees of freedom
Multiple R-squared:  0.821, Adjusted R-squared:  0.781 
F-statistic: 20.6 on 2 and 9 DF,  p-value: 0.000438

tab_model(m4, digits=3, show.ci=FALSE,
          show.se=TRUE, collapse.se=TRUE)

	Y
Predictors	Estimates	p
(Intercept)	0.250 (0.640)	0.705
A	1.500 (0.739)	0.073
B	4.500 (0.739)	<0.001
Observations	12
R2 / R2 adjusted	0.821 / 0.781

## 截距含义：A=B=0时，Y的预测期望值（此时不再是简单均值）
mean(data[A==0 & B==0]$Y)  # 简单均值

[1] 1

predict(m4, data.table(A=0, B=0))  # A=B=0时，Y的预测期望值

   1 
0.25 

## 斜率含义：A或B从0变到1，Y的预测期望值变化
mean(data[A==1]$Y) - mean(data[A==0]$Y)

[1] 1.5

mean(data[B==1]$Y) - mean(data[B==0]$Y)

[1] 4.5

【实践4】交互作用模型

• R语言回归模型公式中的交互作用设定（适用于所有类型的回归模型）
  • Y ~ A * B（精简写法）
  • Y ~ 1 + A + B + A:B（展开写法）
    • 1为截距（可省略不写）
    • A和B为固定效应（“条件效应”，并非ANOVA中的“主效应”！！！）⭐️
    • A:B为交互项（interaction term）

## 方差分析（回顾）
MANOVA(data, dv="Y", between=c("A", "B"))

====== ANOVA (Between-Subjects Design) ======

Descriptives:
────────────────────────
 "A" "B"  Mean    S.D. n
────────────────────────
  A0  B0 1.000 (1.000) 3
  A0  B1 4.000 (1.000) 3
  A1  B0 1.000 (1.000) 3
  A1  B1 7.000 (1.000) 3
────────────────────────
Total sample size: N = 12

ANOVA Table:
Dependent variable(s):      Y
Between-subjects factor(s): A, B
Within-subjects factor(s):  –
Covariate(s):               –
─────────────────────────────────────────────────────────────────────
           MS   MSE df1 df2      F     p     η²p [90% CI of η²p]  η²G
─────────────────────────────────────────────────────────────────────
A       6.750 1.000   1   8  6.750  .032 *     .458 [.034, .715] .458
B      60.750 1.000   1   8 60.750 <.001 ***   .884 [.680, .940] .884
A * B   6.750 1.000   1   8  6.750  .032 *     .458 [.034, .715] .458
─────────────────────────────────────────────────────────────────────
MSE = mean square error (the residual variance of the linear model)
η²p = partial eta-squared = SS / (SS + SSE) = F * df1 / (F * df1 + df2)
ω²p = partial omega-squared = (F - 1) * df1 / (F * df1 + df2 + 1)
η²G = generalized eta-squared (see Olejnik & Algina, 2003)
Cohen’s f² = η²p / (1 - η²p)

Levene’s Test for Homogeneity of Variance:
───────────────────────────────────
       Levene’s F df1 df2     p    
───────────────────────────────────
DV: Y       0.000   3   8 1.000    
───────────────────────────────────

## 回归分析
m5 = lm(Y ~ A * B, data)  # 相当于：Y ~ 1 + A + B + A:B
anova(m5)  # 真正的主效应、交互作用

Analysis of Variance Table

Response: Y
          Df Sum Sq Mean Sq F value  Pr(>F)    
A          1    6.7     6.7    6.75   0.032 *  
B          1   60.7    60.7   60.75 5.3e-05 ***
A:B        1    6.8     6.8    6.75   0.032 *  
Residuals  8    8.0     1.0                    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

summary(m5)  # A和B的固定效应（条件效应），并非“主效应”！！！

Call:
lm(formula = Y ~ A * B, data = data)

Residuals:
   Min     1Q Median     3Q    Max 
    -1     -1      0      1      1 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)   
(Intercept)  1.00e+00   5.77e-01    1.73   0.1215   
A           -1.09e-15   8.16e-01    0.00   1.0000   
B            3.00e+00   8.16e-01    3.67   0.0063 **
A:B          3.00e+00   1.15e+00    2.60   0.0317 * 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1 on 8 degrees of freedom
Multiple R-squared:  0.903, Adjusted R-squared:  0.866 
F-statistic: 24.8 on 3 and 8 DF,  p-value: 0.000212

## A和B的固定效应（条件效应），并非“主效应”！！！
tab_model(m5, digits=3, show.ci=FALSE,
          show.se=TRUE, collapse.se=TRUE)

	Y
Predictors	Estimates	p
(Intercept)	1.000 (0.577)	0.122
A	-0.000 (0.816)	1.000
B	3.000 (0.816)	0.006
A × B	3.000 (1.155)	0.032
Observations	12
R2 / R2 adjusted	0.903 / 0.866

【知识点】为什么回归分析交互作用中的斜率并不是“主效应”？

• 回归分析中，如果存在交互作用（A × B），则每个自变量的斜率代表另一个变量取0时，该变量的简单斜率（也可以叫：固定效应、简单效应、条件效应）！
• 在【实践4】模型结果中：
  • A的斜率/条件效应（b = 0.00）：B = 0时，A的简单斜率（1.00 – 1.00 = 0.00）
  • B的斜率/条件效应（b = 3.00）：A = 0时，B的简单斜率（4.00 – 1.00 = 3.00）
• 而真正的主效应结果是：
  • A的主效应/边际效应：无论B，A1与A0的边际均值差异（4.00 – 2.50 = 1.50）
  • B的主效应/边际效应：无论A，B1与B0的边际均值差异（5.50 – 1.00 = 4.50）

表1. 原始数据

原始数据	B0条件	B1条件
A0条件	0, 1, 2	3, 4, 5
A1条件	0, 1, 2	6, 7, 8

表2. A和B的条件均值与边际均值

均值	B0条件均值	B1条件均值	边际均值
A0条件均值	1.00	4.00	2.50
A1条件均值	1.00	7.00	4.00
边际均值	1.00	5.50	3.25

## 注意，这不是A和B的主效应！（这是固定效应、简单效应、条件效应！）
tab_model(m5, digits=3, show.ci=FALSE,
          show.se=TRUE, collapse.se=TRUE)

	Y
Predictors	Estimates	p
(Intercept)	1.000 (0.577)	0.122
A	-0.000 (0.816)	1.000
B	3.000 (0.816)	0.006
A × B	3.000 (1.155)	0.032
Observations	12
R2 / R2 adjusted	0.903 / 0.866

## 对比方差分析“简单效应”结果
EMMEANS(m, "A", by="B")

------ EMMEANS (effect = "A") ------

Joint Tests of "A":
────────────────────────────────────────────────────────
 Effect "B" df1 df2      F     p     η²p [90% CI of η²p]
────────────────────────────────────────────────────────
      A  B0   1   8  0.000 1.000       .000 [.000, .000]
      A  B1   1   8 13.500  .006 **    .628 [.191, .808]
────────────────────────────────────────────────────────
Note. Simple effects of repeated measures with 3 or more levels
are different from the results obtained with SPSS MANOVA syntax.

Univariate Tests of "A":
─────────────────────────────────────────────────────────
           Sum of Squares df Mean Square      F     p    
─────────────────────────────────────────────────────────
B0: "A"             0.000  1       0.000  0.000 1.000    
B1: "A"            13.500  1      13.500 13.500  .006 ** 
Residuals           8.000  8       1.000                 
─────────────────────────────────────────────────────────
Note. Identical to the results obtained with SPSS GLM EMMEANS syntax.

Estimated Marginal Means of "A":
──────────────────────────────────────
 "A" "B" Mean [95% CI of Mean]    S.E.
──────────────────────────────────────
  A0  B0 1.000 [-0.331, 2.331] (0.577)
  A1  B0 1.000 [-0.331, 2.331] (0.577)
  A0  B1 4.000 [ 2.669, 5.331] (0.577)
  A1  B1 7.000 [ 5.669, 8.331] (0.577)
──────────────────────────────────────

Pairwise Comparisons of "A":
─────────────────────────────────────────────────────────────────────────
 Contrast "B" Estimate    S.E. df     t     p     Cohen’s d [95% CI of d]
─────────────────────────────────────────────────────────────────────────
  A1 - A0  B0    0.000 (0.816)  8 0.000 1.000       0.000 [-1.883, 1.883]
  A1 - A0  B1    3.000 (0.816)  8 3.674  .006 **    3.000 [ 1.117, 4.883]
─────────────────────────────────────────────────────────────────────────
Pooled SD for computing Cohen’s d: 1.000
No need to adjust p values.

Disclaimer:
By default, pooled SD is Root Mean Square Error (RMSE).
There is much disagreement on how to compute Cohen’s d.
You are completely responsible for setting `sd.pooled`.
You might also use `effectsize::t_to_d()` to compute d.

【实践5】如何让回归分析交互作用中的斜率具有“主效应”含义？

• 变量中心化（Chap05）：使新变量（B’）取0等价于原变量（B）取均值
  • 中心化消除了“边际效应”与“条件效应”之间的偏差
    • 边际效应（主效应）= 条件效应（固定效应/回归斜率）+ 交互效应 × 调节变量均值
    • \(\text{Main Effect}_X = b_{X_{Marginal}} = b_{X_{Conditional}} + b_{X*M_{Interaction}} \bar{M}\)
    • \(\text{Main Effect}_A = 1.50 = 0.00 + 3.00 \times (0+1)/2\)
    • \(\text{Main Effect}_B = 4.50 = 3.00 + 3.00 \times (0+1)/2\)

## 中心化处理
data$A.c = data$A - mean(data$A)  # 也可以直接用 grand_mean_center()
data$B.c = data$B - mean(data$B)  # 也可以直接用 grand_mean_center()
data

索引(index): <A>, <B>, <B__A>
        A     B     Y   A.c   B.c
    <num> <num> <num> <num> <num>
 1:     0     0     0  -0.5  -0.5
 2:     0     0     1  -0.5  -0.5
 3:     0     0     2  -0.5  -0.5
 4:     0     1     3  -0.5   0.5
 5:     0     1     4  -0.5   0.5
 6:     0     1     5  -0.5   0.5
 7:     1     0     0   0.5  -0.5
 8:     1     0     1   0.5  -0.5
 9:     1     0     2   0.5  -0.5
10:     1     1     6   0.5   0.5
11:     1     1     7   0.5   0.5
12:     1     1     8   0.5   0.5

m6 = lm(Y ~ A.c * B.c, data)
anova(m6)  # 真正的主效应、交互作用

Analysis of Variance Table

Response: Y
          Df Sum Sq Mean Sq F value  Pr(>F)    
A.c        1    6.8     6.8    6.75   0.032 *  
B.c        1   60.7    60.7   60.75 5.3e-05 ***
A.c:B.c    1    6.8     6.8    6.75   0.032 *  
Residuals  8    8.0     1.0                    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

summary(m6)  # A和B中心化后，固定效应（条件效应）等价于“主效应”

Call:
lm(formula = Y ~ A.c * B.c, data = data)

Residuals:
   Min     1Q Median     3Q    Max 
    -1     -1      0      1      1 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)    3.250      0.289   11.26  3.5e-06 ***
A.c            1.500      0.577    2.60    0.032 *  
B.c            4.500      0.577    7.79  5.3e-05 ***
A.c:B.c        3.000      1.155    2.60    0.032 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1 on 8 degrees of freedom
Multiple R-squared:  0.903, Adjusted R-squared:  0.866 
F-statistic: 24.8 on 3 and 8 DF,  p-value: 0.000212

## A和B中心化后，固定效应（条件效应）等价于“主效应”
tab_model(m6, digits=3, show.ci=FALSE,
          show.se=TRUE, collapse.se=TRUE)

	Y
Predictors	Estimates	p
(Intercept)	3.250 (0.289)	<0.001
A c	1.500 (0.577)	0.032
B c	4.500 (0.577)	<0.001
A c × B c	3.000 (1.155)	0.032
Observations	12
R2 / R2 adjusted	0.903 / 0.866

多个模型汇总与表格输出

【实践6】多元回归模型的结果表输出

模型结果表可以输出到：

• R Console控制台
• Microsoft Word文档
• HTML表格

m1 = lm(Y ~ 1, data)
m2 = lm(Y ~ A, data)
m3 = lm(Y ~ B, data)
m4 = lm(Y ~ A + B, data)
m5 = lm(Y ~ A * B, data)
m6 = lm(Y ~ A.c * B.c, data)

model_summary(list(m1, m2, m3, m4, m5, m6))

Model Summary

───────────────────────────────────────────────────────────────────────────────
             (1) Y      (2) Y     (3) Y       (4) Y       (5) Y      (6) Y     
───────────────────────────────────────────────────────────────────────────────
(Intercept)   3.250 **   2.500 *   1.000       0.250       1.000      3.250 ***
             (0.789)    (1.122)   (0.599)     (0.640)     (0.577)    (0.289)   
A                        1.500                 1.500      -0.000               
                        (1.586)               (0.739)     (0.816)              
B                                  4.500 ***   4.500 ***   3.000 **            
                                  (0.847)     (0.739)     (0.816)              
A:B                                                        3.000 *             
                                                          (1.155)              
A.c                                                                   1.500 *  
                                                                     (0.577)   
B.c                                                                   4.500 ***
                                                                     (0.577)   
A.c:B.c                                                               3.000 *  
                                                                     (1.155)   
───────────────────────────────────────────────────────────────────────────────
R^2           0.000      0.082     0.739       0.821       0.903      0.903    
Adj. R^2      0.000     -0.010     0.712       0.781       0.866      0.866    
Num. obs.    12         12        12          12          12         12        
───────────────────────────────────────────────────────────────────────────────
Note. * p < .05, ** p < .01, *** p < .001.

model_summary(list(m1, m2, m3, m4, m5, m6),
              digits=2,
              file="Models.doc")

✔ Table saved to 'D:/RStudio/RCourse/Models.doc'

tab_model(m1, m2, m3, m4, m5, m6,
          digits=2, string.est="b (SE)",
          show.ci=FALSE, show.se=TRUE, collapse.se=TRUE)

	Y		Y		Y		Y		Y		Y
Predictors	b (SE)	p	b (SE)	p	b (SE)	p	b (SE)	p	b (SE)	p	b (SE)	p
(Intercept)	3.25 (0.79)	0.002	2.50 (1.12)	0.050	1.00 (0.60)	0.126	0.25 (0.64)	0.705	1.00 (0.58)	0.122	3.25 (0.29)	<0.001
A			1.50 (1.59)	0.367			1.50 (0.74)	0.073	-0.00 (0.82)	1.000
B					4.50 (0.85)	<0.001	4.50 (0.74)	<0.001	3.00 (0.82)	0.006
A × B									3.00 (1.15)	0.032
A c											1.50 (0.58)	0.032
B c											4.50 (0.58)	<0.001
A c × B c											3.00 (1.15)	0.032
Observations	12		12		12		12		12		12
R2 / R2 adjusted	0.000 / 0.000		0.082 / -0.010		0.739 / 0.712		0.821 / 0.781		0.903 / 0.866		0.903 / 0.866

【实践7】空气质量数据回归分析实例

## 数据准备
str(airquality)

'data.frame':   153 obs. of  6 variables:
 $ Ozone  : int  41 36 12 18 NA 28 23 19 8 NA ...
 $ Solar.R: int  190 118 149 313 NA NA 299 99 19 194 ...
 $ Wind   : num  7.4 8 12.6 11.5 14.3 14.9 8.6 13.8 20.1 8.6 ...
 $ Temp   : int  67 72 74 62 56 66 65 59 61 69 ...
 $ Month  : int  5 5 5 5 5 5 5 5 5 5 ...
 $ Day    : int  1 2 3 4 5 6 7 8 9 10 ...

Corr(airquality)

Pearson's r and 95% confidence intervals:
─────────────────────────────────────────────────
                   r       [95% CI]     p       N
─────────────────────────────────────────────────
Ozone-Solar.R   0.35 [ 0.17,  0.50] <.001 *** 111
Ozone-Wind     -0.60 [-0.71, -0.47] <.001 *** 116
Ozone-Temp      0.70 [ 0.59,  0.78] <.001 *** 116
Ozone-Month     0.16 [-0.02,  0.34]  .078 .   116
Ozone-Day      -0.01 [-0.20,  0.17]  .888     116
Solar.R-Wind   -0.06 [-0.22,  0.11]  .496     146
Solar.R-Temp    0.28 [ 0.12,  0.42] <.001 *** 146
Solar.R-Month  -0.08 [-0.23,  0.09]  .366     146
Solar.R-Day    -0.15 [-0.31,  0.01]  .070 .   146
Wind-Temp      -0.46 [-0.57, -0.32] <.001 *** 153
Wind-Month     -0.18 [-0.33, -0.02]  .027 *   153
Wind-Day        0.03 [-0.13,  0.19]  .739     153
Temp-Month      0.42 [ 0.28,  0.54] <.001 *** 153
Temp-Day       -0.13 [-0.28,  0.03]  .108     153
Month-Day      -0.01 [-0.17,  0.15]  .922     153
─────────────────────────────────────────────────

## 回归分析
m1 = lm(Ozone ~ Solar.R, data=airquality)
m2 = lm(Ozone ~ Solar.R + Temp, data=airquality)
m3 = lm(Ozone ~ Solar.R + Temp + Wind, data=airquality)
m4 = lm(Ozone ~ Solar.R + Temp + Wind + Month, data=airquality)
model = lm(Ozone ~ ., data=airquality)  # . 表示纳入剩余所有变量

思考：

• 截距（b = –64.116, p = .007）的含义是？
• 气温斜率（b = 1.896, β = .543, p < .001）的含义是？
• t值与b、SE的关系是？（请验算一下）

\[ t = \frac{b}{SE} \]

## bruceR系列函数
GLM_summary(model)

General Linear Model (OLS Regression)

Model Fit:
F(5, 105) = 34.99, p = 7e-21 ***
R² = 0.62494 (Adjusted R² = 0.60708)

Unstandardized Coefficients:
Outcome Variable: Ozone
N = 111 (42 missing cases deleted)
────────────────────────────────────────────────────────────────────────
                   b     S.E.      t     p           [95% CI of b]   VIF
────────────────────────────────────────────────────────────────────────
(Intercept)  -64.116 (23.482) -2.730  .007 **  [-110.678, -17.555]      
Solar.R        0.050 ( 0.023)  2.147  .034 *   [   0.004,   0.097] 1.152
Wind          -3.318 ( 0.645) -5.149 <.001 *** [  -4.596,  -2.041] 1.329
Temp           1.896 ( 0.274)  6.922 <.001 *** [   1.353,   2.439] 1.722
Month         -3.040 ( 1.513) -2.009  .047 *   [  -6.041,  -0.039] 1.257
Day            0.274 ( 0.230)  1.192  .236     [  -0.182,   0.729] 1.011
────────────────────────────────────────────────────────────────────────

Standardized Coefficients (β):
Outcome Variable: Ozone
N = 111 (42 missing cases deleted)
────────────────────────────────────────────────────────────────────────────
              β    S.E.      t     p        [95% CI of β] r(partial) r(part)
────────────────────────────────────────────────────────────────────────────
Solar.R   0.138 (0.064)  2.147  .034 *   [ 0.011,  0.265]      0.205   0.128
Wind     -0.355 (0.069) -5.149 <.001 *** [-0.491, -0.218]     -0.449  -0.308
Temp      0.543 (0.078)  6.922 <.001 *** [ 0.387,  0.698]      0.560   0.414
Month    -0.135 (0.067) -2.009  .047 *   [-0.267, -0.002]     -0.192  -0.120
Day       0.072 (0.060)  1.192  .236     [-0.047,  0.191]      0.116   0.071
────────────────────────────────────────────────────────────────────────────

model_summary(model)

Model Summary

────────────────────────
             (1) Ozone  
────────────────────────
(Intercept)  -64.116 ** 
             (23.482)   
Solar.R        0.050 *  
              (0.023)   
Wind          -3.318 ***
              (0.645)   
Temp           1.896 ***
              (0.274)   
Month         -3.040 *  
              (1.513)   
Day            0.274    
              (0.230)   
────────────────────────
R^2            0.625    
Adj. R^2       0.607    
Num. obs.    111        
────────────────────────
Note. * p < .05, ** p < .01, *** p < .001.

# Check for Multicollinearity

Low Correlation

    Term  VIF       VIF 95% CI adj. VIF Tolerance Tolerance 95% CI
 Solar.R 1.15 [1.03,     1.66]     1.07      0.87     [0.60, 0.97]
    Wind 1.33 [1.14,     1.78]     1.15      0.75     [0.56, 0.88]
    Temp 1.72 [1.40,     2.29]     1.31      0.58     [0.44, 0.71]
   Month 1.26 [1.09,     1.71]     1.12      0.80     [0.59, 0.91]
     Day 1.01 [1.00, 83489.99]     1.01      0.99     [0.00, 1.00]

model_summary(model, std=TRUE)  # 标准化回归系数

Model Summary

──────────────────────
           (1) Ozone  
──────────────────────
Solar.R       .138 *  
             (.064)   
Wind         -.355 ***
             (.069)   
Temp          .543 ***
             (.078)   
Month        -.135 *  
             (.067)   
Day           .072    
             (.060)   
──────────────────────
R^2           .625    
Adj. R^2      .607    
Num. obs.  111        
──────────────────────
Note. * p < .05, ** p < .01, *** p < .001.

# Check for Multicollinearity

Low Correlation

    Term  VIF       VIF 95% CI adj. VIF Tolerance Tolerance 95% CI
 Solar.R 1.15 [1.03,     1.66]     1.07      0.87     [0.60, 0.97]
    Wind 1.33 [1.14,     1.78]     1.15      0.75     [0.56, 0.88]
    Temp 1.72 [1.40,     2.29]     1.31      0.58     [0.44, 0.71]
   Month 1.26 [1.09,     1.71]     1.12      0.80     [0.59, 0.91]
     Day 1.01 [1.00, 83489.99]     1.01      0.99     [0.00, 1.00]

model_summary(list(m1, m2, m3, m4, model), file="Models.doc")

✔ Table saved to 'D:/RStudio/RCourse/Models.doc'

↑ APA格式表格，可直接复制到论文中

## sjPlot系列函数
tab_model(model)

	Ozone
Predictors	Estimates	CI	p
(Intercept)	-64.12	-110.68 – -17.55	0.007
Solar R	0.05	0.00 – 0.10	0.034
Wind	-3.32	-4.60 – -2.04	<0.001
Temp	1.90	1.35 – 2.44	<0.001
Month	-3.04	-6.04 – -0.04	0.047
Day	0.27	-0.18 – 0.73	0.236
Observations	111
R2 / R2 adjusted	0.625 / 0.607

tab_model(m1, m2, m3, m4, model,
          digits=3, string.est="b (SE)",
          show.ci=FALSE, show.se=TRUE, collapse.se=TRUE)

	Ozone		Ozone		Ozone		Ozone		Ozone
Predictors	b (SE)	p	b (SE)	p	b (SE)	p	b (SE)	p	b (SE)	p
(Intercept)	18.599 (6.748)	0.007	-145.703 (18.447)	<0.001	-64.342 (23.055)	0.006	-58.054 (22.971)	0.013	-64.116 (23.482)	0.007
Solar R	0.127 (0.033)	<0.001	0.057 (0.026)	0.028	0.060 (0.023)	0.011	0.050 (0.023)	0.037	0.050 (0.023)	0.034
Temp			2.278 (0.246)	<0.001	1.652 (0.254)	<0.001	1.871 (0.274)	<0.001	1.896 (0.274)	<0.001
Wind					-3.334 (0.654)	<0.001	-3.317 (0.646)	<0.001	-3.318 (0.645)	<0.001
Month							-2.992 (1.516)	0.051	-3.040 (1.513)	0.047
Day									0.274 (0.230)	0.236
Observations	111		111		111		111		111
R2 / R2 adjusted	0.121 / 0.113		0.510 / 0.501		0.606 / 0.595		0.620 / 0.606		0.625 / 0.607

简单斜率检验与多重比较

【实践8】基于PROCESS()函数的调节效应与简单斜率分析

• bruceR包PROCESS()函数：中介与调节效应分析多功能“一站式”函数
  • PROCESS()函数帮助文档
  • PROCESS()函数轻松实现中介效应和调节效应分析 - 知乎文章
• 调节效应 = 交互作用
  • 调节效应本质上就是交互作用，只是明确了交互项中的哪个变量是调节变量（moderator）

## 数据准备
data = airquality
data$Month = as.factor(data$Month)
data$Temp.C = (data$Temp - 32) / 1.8  # 摄氏度 = (华氏度 - 32) / 1.8
str(data)

'data.frame':   153 obs. of  7 variables:
 $ Ozone  : int  41 36 12 18 NA 28 23 19 8 NA ...
 $ Solar.R: int  190 118 149 313 NA NA 299 99 19 194 ...
 $ Wind   : num  7.4 8 12.6 11.5 14.3 14.9 8.6 13.8 20.1 8.6 ...
 $ Temp   : int  67 72 74 62 56 66 65 59 61 69 ...
 $ Month  : Factor w/ 5 levels "5","6","7","8",..: 1 1 1 1 1 1 1 1 1 1 ...
 $ Day    : int  1 2 3 4 5 6 7 8 9 10 ...
 $ Temp.C : num  19.4 22.2 23.3 16.7 13.3 ...

## 描述统计
Freq(data$Month)

Frequency Statistics:
──────────
    N    %
──────────
5  31 20.3
6  30 19.6
7  31 20.3
8  31 20.3
9  30 19.6
──────────
Total N = 153

Describe(data)

NOTE: `Month` transformed to numeric.

Descriptive Statistics:
──────────────────────────────────────────────────────────────────────
           N (NA)   Mean    SD | Median   Min    Max Skewness Kurtosis
──────────────────────────────────────────────────────────────────────
Ozone    116   37  42.13 32.99 |  31.50  1.00 168.00     1.21     1.11
Solar.R  146    7 185.93 90.06 | 205.00  7.00 334.00    -0.42    -1.00
Wind     153        9.96  3.52 |   9.70  1.70  20.70     0.34     0.03
Temp     153       77.88  9.47 |  79.00 56.00  97.00    -0.37    -0.46
Month*   153        2.99  1.42 |   3.00  1.00   5.00    -0.00    -1.32
Day      153       15.80  8.86 |  16.00  1.00  31.00     0.00    -1.22
Temp.C   153       25.49  5.26 |  26.11 13.33  36.11    -0.37    -0.46
──────────────────────────────────────────────────────────────────────

## PROCESS()调节效应分析
PROCESS(data, y="Ozone", x="Solar.R", mods="Month")

****************** PART 1. Regression Model Summary ******************

PROCESS Model ID : 1
Model Type : Simple Moderation
-    Outcome (Y) : Ozone
-  Predictor (X) : Solar.R
-  Mediators (M) : -
- Moderators (W) : Month
- Covariates (C) : -
-   HLM Clusters : -

All numeric predictors have been grand-mean centered.
(For details, please see the help page of PROCESS.)

Formula of Outcome:
-    Ozone ~ Solar.R*Month

CAUTION:
  Fixed effect (coef.) of a predictor involved in an interaction
  denotes its "simple effect/slope" at the other predictor = 0.
  Only when all predictors in an interaction are mean-centered
  can the fixed effect be interpreted as "main effect"!
  
Model Summary

────────────────────────────────────────
                (1) Ozone    (2) Ozone  
────────────────────────────────────────
(Intercept)      42.099 ***   24.254 ***
                 (2.974)      (5.617)   
Solar.R           0.127 ***    0.047    
                 (0.033)      (0.048)   
Month6                         5.269    
                             (10.753)   
Month7                        29.525 ***
                              (8.082)   
Month8                        39.033 ***
                              (8.075)   
Month9                         8.093    
                              (7.668)   
Solar.R:Month6                 0.090    
                              (0.113)   
Solar.R:Month7                 0.122    
                              (0.084)   
Solar.R:Month8                 0.234 ** 
                              (0.089)   
Solar.R:Month9                 0.007    
                              (0.081)   
────────────────────────────────────────
R^2               0.121        0.372    
Adj. R^2          0.113        0.317    
Num. obs.       111          111        
────────────────────────────────────────
Note. * p < .05, ** p < .01, *** p < .001.

************ PART 2. Mediation/Moderation Effect Estimate ************

Package Use : ‘interactions’ (v1.2.0)
Effect Type : Simple Moderation (Model 1)
Sample Size : 111 (42 missing observations deleted)
Random Seed : -
Simulations : -

Interaction Effect on "Ozone" (Y)
───────────────────────────────────────
                    F df1 df2     p    
───────────────────────────────────────
Solar.R * Month  2.12   4 101  .083 .  
───────────────────────────────────────

Simple Slopes: "Solar.R" (X) ==> "Ozone" (Y)
───────────────────────────────────────────────────────
 "Month" Effect    S.E.     t     p            [95% CI]
───────────────────────────────────────────────────────
 5        0.047 (0.048) 0.969  .335     [-0.049, 0.143]
 6        0.137 (0.102) 1.344  .182     [-0.065, 0.338]
 7        0.169 (0.068) 2.469  .015 *   [ 0.033, 0.304]
 8        0.281 (0.074) 3.772 <.001 *** [ 0.133, 0.428]
 9        0.054 (0.065) 0.838  .404     [-0.074, 0.182]
───────────────────────────────────────────────────────

PROCESS(data, y="Ozone", x="Solar.R", mods="Temp.C")

****************** PART 1. Regression Model Summary ******************

PROCESS Model ID : 1
Model Type : Simple Moderation
-    Outcome (Y) : Ozone
-  Predictor (X) : Solar.R
-  Mediators (M) : -
- Moderators (W) : Temp.C
- Covariates (C) : -
-   HLM Clusters : -

All numeric predictors have been grand-mean centered.
(For details, please see the help page of PROCESS.)

Formula of Outcome:
-    Ozone ~ Solar.R*Temp.C

CAUTION:
  Fixed effect (coef.) of a predictor involved in an interaction
  denotes its "simple effect/slope" at the other predictor = 0.
  Only when all predictors in an interaction are mean-centered
  can the fixed effect be interpreted as "main effect"!
  
Model Summary

────────────────────────────────────────
                (1) Ozone    (2) Ozone  
────────────────────────────────────────
(Intercept)      42.099 ***   40.004 ***
                 (2.974)      (2.247)   
Solar.R           0.127 ***    0.096 ***
                 (0.033)      (0.028)   
Temp.C                         4.122 ***
                              (0.426)   
Solar.R:Temp.C                 0.015 ** 
                              (0.005)   
────────────────────────────────────────
R^2               0.121        0.551    
Adj. R^2          0.113        0.539    
Num. obs.       111          111        
────────────────────────────────────────
Note. * p < .05, ** p < .01, *** p < .001.

************ PART 2. Mediation/Moderation Effect Estimate ************

Package Use : ‘interactions’ (v1.2.0)
Effect Type : Simple Moderation (Model 1)
Sample Size : 111 (42 missing observations deleted)
Random Seed : -
Simulations : -

Interaction Effect on "Ozone" (Y)
────────────────────────────────────────
                     F df1 df2     p    
────────────────────────────────────────
Solar.R * Temp.C  9.77   1 107  .002 ** 
────────────────────────────────────────

Simple Slopes: "Solar.R" (X) ==> "Ozone" (Y)
─────────────────────────────────────────────────────────────
 "Temp.C"      Effect    S.E.     t     p            [95% CI]
─────────────────────────────────────────────────────────────
 20.146 (- SD)  0.017 (0.028) 0.613  .541     [-0.038, 0.072]
 25.440 (Mean)  0.096 (0.028) 3.466 <.001 *** [ 0.041, 0.151]
 30.735 (+ SD)  0.175 (0.045) 3.880 <.001 *** [ 0.085, 0.264]
─────────────────────────────────────────────────────────────

PROCESS(data, y="Ozone", x="Solar.R", mods="Temp.C",
        mod1.val=c(0, 5, 10, 15, 20, 25, 30, 35, 40))

****************** PART 1. Regression Model Summary ******************

PROCESS Model ID : 1
Model Type : Simple Moderation
-    Outcome (Y) : Ozone
-  Predictor (X) : Solar.R
-  Mediators (M) : -
- Moderators (W) : Temp.C
- Covariates (C) : -
-   HLM Clusters : -

All numeric predictors have been grand-mean centered.
(For details, please see the help page of PROCESS.)

Formula of Outcome:
-    Ozone ~ Solar.R*Temp.C

CAUTION:
  Fixed effect (coef.) of a predictor involved in an interaction
  denotes its "simple effect/slope" at the other predictor = 0.
  Only when all predictors in an interaction are mean-centered
  can the fixed effect be interpreted as "main effect"!
  
Model Summary

────────────────────────────────────────
                (1) Ozone    (2) Ozone  
────────────────────────────────────────
(Intercept)      42.099 ***   40.004 ***
                 (2.974)      (2.247)   
Solar.R           0.127 ***    0.096 ***
                 (0.033)      (0.028)   
Temp.C                         4.122 ***
                              (0.426)   
Solar.R:Temp.C                 0.015 ** 
                              (0.005)   
────────────────────────────────────────
R^2               0.121        0.551    
Adj. R^2          0.113        0.539    
Num. obs.       111          111        
────────────────────────────────────────
Note. * p < .05, ** p < .01, *** p < .001.

************ PART 2. Mediation/Moderation Effect Estimate ************

Package Use : ‘interactions’ (v1.2.0)
Effect Type : Simple Moderation (Model 1)
Sample Size : 111 (42 missing observations deleted)
Random Seed : -
Simulations : -

Interaction Effect on "Ozone" (Y)
────────────────────────────────────────
                     F df1 df2     p    
────────────────────────────────────────
Solar.R * Temp.C  9.77   1 107  .002 ** 
────────────────────────────────────────

Simple Slopes: "Solar.R" (X) ==> "Ozone" (Y)
──────────────────────────────────────────────────────────
 "Temp.C" Effect    S.E.      t     p             [95% CI]
──────────────────────────────────────────────────────────
 0.000    -0.283 (0.112) -2.536  .013 *   [-0.504, -0.062]
 5.000    -0.209 (0.089) -2.356  .020 *   [-0.384, -0.033]
 10.000   -0.134 (0.066) -2.032  .045 *   [-0.265, -0.003]
 15.000   -0.060 (0.045) -1.330  .186     [-0.148,  0.029]
 20.000    0.015 (0.028)  0.529  .598     [-0.041,  0.071]
 25.000    0.089 (0.027)  3.335  .001 **  [ 0.036,  0.143]
 30.000    0.164 (0.042)  3.886 <.001 *** [ 0.080,  0.247]
 35.000    0.238 (0.063)  3.781 <.001 *** [ 0.113,  0.363]
 40.000    0.313 (0.085)  3.660 <.001 *** [ 0.143,  0.482]
──────────────────────────────────────────────────────────

PROCESS(data, y="Ozone", x="Solar.R", mods="Temp.C",
        covs="Month")  # Month作为协变量（控制变量）

****************** PART 1. Regression Model Summary ******************

PROCESS Model ID : 1
Model Type : Simple Moderation
-    Outcome (Y) : Ozone
-  Predictor (X) : Solar.R
-  Mediators (M) : -
- Moderators (W) : Temp.C
- Covariates (C) : Month
-   HLM Clusters : -

All numeric predictors have been grand-mean centered.
(For details, please see the help page of PROCESS.)

Formula of Outcome:
-    Ozone ~ Month + Solar.R*Temp.C

CAUTION:
  Fixed effect (coef.) of a predictor involved in an interaction
  denotes its "simple effect/slope" at the other predictor = 0.
  Only when all predictors in an interaction are mean-centered
  can the fixed effect be interpreted as "main effect"!
  
Model Summary

────────────────────────────────────────
                (1) Ozone    (2) Ozone  
────────────────────────────────────────
(Intercept)      24.442 ***   49.873 ***
                 (5.735)      (5.886)   
Month6            5.069      -24.915 ** 
                (10.980)      (9.370)   
Month7           31.042 ***  -11.393    
                 (8.018)      (8.254)   
Month8           36.903 ***   -5.166    
                 (8.201)      (8.607)   
Month9            8.912      -15.473 *  
                 (7.763)      (7.080)   
Solar.R           0.115 ***    0.088 ** 
                 (0.030)      (0.028)   
Temp.C                         4.553 ***
                              (0.588)   
Solar.R:Temp.C                 0.014 ** 
                              (0.005)   
────────────────────────────────────────
R^2               0.320        0.595    
Adj. R^2          0.287        0.567    
Num. obs.       111          111        
────────────────────────────────────────
Note. * p < .05, ** p < .01, *** p < .001.

************ PART 2. Mediation/Moderation Effect Estimate ************

Package Use : ‘interactions’ (v1.2.0)
Effect Type : Simple Moderation (Model 1)
Sample Size : 111 (42 missing observations deleted)
Random Seed : -
Simulations : -

Interaction Effect on "Ozone" (Y)
────────────────────────────────────────
                     F df1 df2     p    
────────────────────────────────────────
Solar.R * Temp.C  9.24   1 103  .003 ** 
────────────────────────────────────────

Simple Slopes: "Solar.R" (X) ==> "Ozone" (Y)
─────────────────────────────────────────────────────────────
 "Temp.C"      Effect    S.E.     t     p            [95% CI]
─────────────────────────────────────────────────────────────
 20.146 (- SD)  0.011 (0.028) 0.416  .679     [-0.043, 0.066]
 25.440 (Mean)  0.088 (0.028) 3.090  .003 **  [ 0.031, 0.144]
 30.735 (+ SD)  0.164 (0.046) 3.566 <.001 *** [ 0.073, 0.255]
─────────────────────────────────────────────────────────────

P = PROCESS(data, y="Ozone", x="Solar.R", mods="Temp.C",
            covs=c("Month", "Day", "Wind"))  # 控制其他所有变量

****************** PART 1. Regression Model Summary ******************

PROCESS Model ID : 1
Model Type : Simple Moderation
-    Outcome (Y) : Ozone
-  Predictor (X) : Solar.R
-  Mediators (M) : -
- Moderators (W) : Temp.C
- Covariates (C) : Month, Day, Wind
-   HLM Clusters : -

All numeric predictors have been grand-mean centered.
(For details, please see the help page of PROCESS.)

Formula of Outcome:
-    Ozone ~ Month + Day + Wind + Solar.R*Temp.C

CAUTION:
  Fixed effect (coef.) of a predictor involved in an interaction
  denotes its "simple effect/slope" at the other predictor = 0.
  Only when all predictors in an interaction are mean-centered
  can the fixed effect be interpreted as "main effect"!
  
Model Summary

────────────────────────────────────────
                (1) Ozone    (2) Ozone  
────────────────────────────────────────
(Intercept)      31.644 ***   47.479 ***
                 (4.927)      (5.394)   
Month6            8.331      -14.987    
                 (9.237)      (8.722)   
Month7           17.823 *     -9.616    
                 (7.005)      (7.531)   
Month8           24.416 ***   -2.389    
                 (7.130)      (7.905)   
Month9            2.088      -13.028 *  
                 (6.594)      (6.444)   
Day               0.055        0.129    
                 (0.261)      (0.223)   
Wind             -4.640 ***   -3.117 ***
                 (0.683)      (0.631)   
Solar.R           0.097 ***    0.093 ***
                 (0.025)      (0.026)   
Temp.C                         3.388 ***
                              (0.592)   
Solar.R:Temp.C                 0.014 ** 
                              (0.004)   
────────────────────────────────────────
R^2               0.530        0.675    
Adj. R^2          0.498        0.646    
Num. obs.       111          111        
────────────────────────────────────────
Note. * p < .05, ** p < .01, *** p < .001.

************ PART 2. Mediation/Moderation Effect Estimate ************

Package Use : ‘interactions’ (v1.2.0)
Effect Type : Simple Moderation (Model 1)
Sample Size : 111 (42 missing observations deleted)
Random Seed : -
Simulations : -

Interaction Effect on "Ozone" (Y)
─────────────────────────────────────────
                      F df1 df2     p    
─────────────────────────────────────────
Solar.R * Temp.C  10.60   1 101  .002 ** 
─────────────────────────────────────────

Simple Slopes: "Solar.R" (X) ==> "Ozone" (Y)
─────────────────────────────────────────────────────────────
 "Temp.C"      Effect    S.E.     t     p            [95% CI]
─────────────────────────────────────────────────────────────
 20.146 (- SD)  0.018 (0.025) 0.717  .475     [-0.032, 0.067]
 25.440 (Mean)  0.093 (0.026) 3.593 <.001 *** [ 0.042, 0.144]
 30.735 (+ SD)  0.167 (0.042) 3.990 <.001 *** [ 0.084, 0.251]
─────────────────────────────────────────────────────────────

P$results[[1]]$jn[[1]]$plot

模型诊断

【知识点】回归模型的基本假定

回归模型的基本假定（basic assumptions）

• 线性关系（linearity）
• 残差独立（independence）
• 残差正态分布（normality）
• 残差方差齐性（homoscedasticity）
• 不存在多重共线性（no multicollinearity）

【实践9】检查回归模型是否满足基本假定

• 使用performance包的check_*()系列函数（bruceR包已加载performance包）

model = lm(Ozone ~ ., data=airquality)

check_model(model)  # 一次性检查所有假定

check_outliers(model)  # 极端值问题

OK: No outliers detected.
- Based on the following method and threshold: cook (0.9).
- For variable: (Whole model)

check_normality(model)  # 非正态问题

Warning: Non-normality of residuals detected (p < .001).

check_heteroscedasticity(model)  # 异方差问题

Warning: Heteroscedasticity (non-constant error variance) detected (p = 0.008).

check_autocorrelation(model)  # 自相关问题

OK: Residuals appear to be independent and not autocorrelated (p = 0.860).

check_collinearity(model)  # 共线性问题

# Check for Multicollinearity

Low Correlation

    Term  VIF       VIF 95% CI adj. VIF Tolerance Tolerance 95% CI
 Solar.R 1.15 [1.03,     1.66]     1.07      0.87     [0.60, 0.97]
    Wind 1.33 [1.14,     1.78]     1.15      0.75     [0.56, 0.88]
    Temp 1.72 [1.40,     2.29]     1.31      0.58     [0.44, 0.71]
   Month 1.26 [1.09,     1.71]     1.12      0.80     [0.59, 0.91]
     Day 1.01 [1.00, 83489.99]     1.01      0.99     [0.00, 1.00]

* 拓展：其他回归模型

【探索发现】各类回归模型使用的R函数

• 一般线性模型
  • lm()函数：OLS回归（连续因变量）
• 广义线性模型
  • glm()函数
    • Logistic回归（二分类因变量）：glm(..., data, family=binomial)
    • Poisson回归（计数因变量）：glm(..., data, family=poisson)

• 多分类因变量回归模型
  • nnet::multinom()函数：Multinomial Log-Linear回归（多分类因变量）
• 线性混合模型（多水平模型）
  • nlme::lme()函数
  • lme4::lmer()、lmerTest::lmer()函数
• 广义线性混合模型
  • nlme::nlme()函数
  • lme4::glmer()函数

（多水平模型/线性混合模型的原理与方法，将在研究生阶段系统学习）

【阶段作业②】统计分析综合

作业要求：

• 基于期末自选数据和所有代码积累，综合运用第7章（基础统计）、第8章（方差分析）、第9章（回归分析）所学的统计方法和代码，自由分析你的期末自选数据，体现对R语言统计分析代码的掌握和迁移应用（7/8/9每章至少有一种方法体现在作业中）
• 使用R Markdown完成，对关键代码及结果要有注释说明

平台提交：

• 运行得到的HTML网页
  • 提交文件命名格式：学号-姓名-R阶段作业2.html