I am trying to get confidence intervals around the predicted summary values after a regression and cannot figure out how to do that in python. I know how to get a prediction and CI for every subject but not for the summary measure. Using the following example my question is: "What is the average gre level for someone with a gpa of 3:"
import pandas as pd
import statsmodels.formula.api as sm
import numpy as np
df=pd.read_stata('https://stats.idre.ucla.edu/stat/stata/dae/binary.dta', convert_categoricals=False)
regresult = sm.ols(formula='gre~gpa', data=df).fit()
pred=regresult.get_prediction(df.assign(gpa=3))
predtable=pred.summary_frame()
print(predtable)
print(np.mean(predtable['mean']))
I know how to get the summary measure [i.e., np.mean(predtable['mean']) ], but not how to get the CIs.
Basically I want to replicate the corresponding Stata output of margins but do not know how:
use http://stats.idre.ucla.edu/stat/stata/dae/binary.dta
regress gre gpa
margins, at(gpa=3)
. margins, at(gpa=3)
RESULT:
Adjusted predictions Number of obs = 400
Model VCE : OLS
Expression : Linear prediction, predict()
at : gpa = 3
------------------------------------------------------------------------------
| Delta-method
| Margin Std. Err. t P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
_cons | 542.2223 7.648634 70.89 0.000 527.1855 557.2591
------------------------------------------------------------------------------
How to get the CI of 527.1855 - 557.2591 in Python?
Best
UPDATE: As Pearly Spencer points out, in the aforementioned example the respective values can be seen in predtable. However, this no longer works when a categorical variable is included:
regresult = sm.ols(formula='gre~gpa+C(rank)', data=df).fit()
pred=regresult.get_prediction(df.assign(gpa=3))
predtable=pred.summary_frame()
print(predtable)
print(np.mean(predtable['mean']))
The Stata commands margins calculates the predicted probability for someone with gpa=3 while setting rank at its mean.
. regress gre gpa i.rank
. margins, at(gpa=3)
Predictive margins Number of obs = 400
Model VCE : OLS
Expression : Linear prediction, predict()
at : gpa = 3
------------------------------------------------------------------------------
| Delta-method
| Margin Std. Err. t P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
_cons | 541.9911 7.640092 70.94 0.000 526.9707 557.0114
------------------------------------------------------------------------------
This is equal to: margins, at(gpa=3) atmeans
In Python print(np.mean(predtable['mean'])) gives the correct value of 541.9911. However, I do not know how to calculate the CIs. Perhaps it would work with using pred=regresult.get_prediction(df.assign(gpa=3)) and somehow include the mean of rank after gpa=3 but I could not get it to work because only rank=1, rank=2, rank=3 and rank=4 are allowed since rank is used as a categorical var in the equation.