使用 python 代码计算 R 中 lm() 函数的输出 "se.fit" 和 "resdiual.scale"

Question

我在 R 中有一个代码片段，其中包含一个简单的 lm() 函数，其中包含一个因变量和一个自变量，如下所示。

X = ([149876.9876, 157853.421, 147822.3803, 147904.6639, 152625.6781, 147229.8083, 181202.081, 164499.6566, 171461.6586, 164309.3919])
Y = ([26212109.07, 28376408.76, 30559566.77, 26765176.65, 28206749.66, 27560521.33, 32713878.83, 31263763.7, 30812063.54, 30225631.6])

lmfit <- lm(formula = Data_df$Y ~ Data_df$X, data=Data_df)
lmpred <- predict(lmfit, newdata=Data_df, se.fit=TRUE, interval = "prediction")
print(lmpred) #prints out fit, se.fit, df, residual.scale

以上代码的输出是3个向量 1.) 适合 2.) se.fit 3.) 自由度 4.) residual.scale

请帮我想办法在python中计算se.fit和residual.scale。 我正在使用 statsmodels.ols 来做线性回归模型。下面是我用来构建线性回归的 python 代码。

import pandas as pd
import statsmodels.formula.api as smf
import numpy as np
ols_result = smf.ols(formula='Y ~ X', data=DATA_X_Y_OLS).fit()
ols_result.predict(data_x_values)

R输出

$fit
        fit      lwr      upr
1  27594475 23262089 31926862
2  28768803 24486082 33051524
3  27291987 22943619 31640354
4  27304101 22956398 31651804
5  27999150 23686118 32312183
6  27204745 22851531 31557960
7  32206302 27951767 36460836
8  29747293 25490577 34004009
9  30772271 26527501 35017042
10 29719281 25462018 33976544


$se.fit
        1         2         3         4         5         6         7         8         9        10        
 578003.4  483363.7  605520.6  604399.0  542961.1  613642.7  420890.0  426036.9  397072.7  427318.3  

$df
[1] 24

$residual.scale
[1] 2017981

Answer 1

为了找到拟合，se.fit、df、residual.scale 是 R 中 lm() 函数的输出。下面是计算上述 4 个值的 python 代码

import statsmodels.formula.api as smf
import numpy as np
ols_result = smf.ols(formula='Y ~ X', data=DATA).fit()
fit = ols_result.predict(X_new) //predicted values ie, fit from lm()

covariance_matrix= ols_result.cov_params()
x = DATA['x'].values
xO = pd.DataFrame({"Constant":np.ones(len(x))}).join(pd.DataFrame(x)).values
x1 = np.dot(xO, COVARIANCE_MATRIX)
se_fit = np.sqrt(np.sum(x1 * xO,axis = 1)) //Standard error of the fitted values ie, se.fit in lm()
df = ols_result.df_resid //Degree of freedom ie, df in lm()

residual_scale = round(np.sqrt(np.dot(np.transpose(x), x)/df)) //Residual SD ie, Residual standard deviation

使用 python 代码计算 R 中 lm() 函数的输出 "se.fit" 和 "resdiual.scale"

To calculate the "se.fit" and "resdiual.scale" that are the outputs from lm() function in R using python code

python

r

linear-regression

lm