import statsmodels.formula.api as ols
X_with_bias = np.append(arr = np.ones((len(X), 1)).astype(int), values = X, axis = 1)
result = ols.OLS(endog = np.array(y).flatten(), exog = X_with_bias).fit()
print(result.summary())

from statsmodels.stats.outliers_influence import variance_inflation_factor
vif = pd.DataFrame()
vif["VIF Factor"] = [variance_inflation_factor(X.values, i) for i in range(X.shape[1])]
vif["features"] = X.columns

from statsmodels.formula.api import ols
regression = ols("data ~ x", data=dict(data=y, x=x)).fit()
test = regression.outlier_test()
outliers = ((x[i],y[i]) for i,t in enumerate(test.icol(2)) if t < 0.5)
print ('Outliers: ', list(outliers))

from statsmodels.compat import lzip
import statsmodels.stats.api as sms
name = ['F statistic', 'p-value']
test = sms.het_goldfeldquandt(result.resid, result.model.exog)
lzip(name, test)

name = ['t value', 'p value']
test = sms.linear_harvey_collier(results)
lzip(name, test)