05 · Coming from R¶

The same analysis in R's performance and in checkmodel, side by side. R code is shown as text; the Python cells run.

# R
library(performance)
m <- lm(mpg ~ wt + hp, data = mtcars)
check_model(m)

In [1]:

import sys

DATA = "../tests/data"
if "google.colab" in sys.modules:  # on Colab: install checkmodel, read the data from GitHub
    %pip install -q checkmodel
    DATA = "https://raw.githubusercontent.com/pablo-legeren/checkmodel/main/tests/data"

import sys DATA = "../tests/data" if "google.colab" in sys.modules: # on Colab: install checkmodel, read the data from GitHub %pip install -q checkmodel DATA = "https://raw.githubusercontent.com/pablo-legeren/checkmodel/main/tests/data"

In [2]:

import warnings

import pandas as pd
import statsmodels.api as sm
import statsmodels.formula.api as smf

warnings.filterwarnings("ignore", category=FutureWarning)
from checkmodel import check_model

mtcars = pd.read_csv(f"{DATA}/mtcars.csv", index_col=0)
m = smf.ols("mpg ~ wt + hp", data=mtcars).fit()
check_model(m).plot();

import warnings import pandas as pd import statsmodels.api as sm import statsmodels.formula.api as smf warnings.filterwarnings("ignore", category=FutureWarning) from checkmodel import check_model mtcars = pd.read_csv(f"{DATA}/mtcars.csv", index_col=0) m = smf.ols("mpg ~ wt + hp", data=mtcars).fit() check_model(m).plot();

check_normality(m)
#> Warning: Non-normality of residuals detected (p = 0.021).
check_heteroscedasticity(m)
#> OK: Error variance appears to be homoscedastic (p = 0.402).

In [3]:

from checkmodel import check_heteroscedasticity, check_normality

print(check_normality(m))
print(check_heteroscedasticity(m))

from checkmodel import check_heteroscedasticity, check_normality print(check_normality(m)) print(check_heteroscedasticity(m))

Warning: Non-normality of residuals detected (p = 0.021).
OK: Error variance appears to be homoscedastic (p = 0.402).

check_collinearity(lm(mpg ~ wt + hp + factor(cyl), data = mtcars))

factor(cyl) is C(cyl) in a Python formula. Terms are listed in formula order, as in R:

In [4]:

from checkmodel import check_collinearity

check_collinearity(smf.ols("mpg ~ wt + hp + C(cyl)", data=mtcars).fit()).to_dataframe().round(3)

from checkmodel import check_collinearity check_collinearity(smf.ols("mpg ~ wt + hp + C(cyl)", data=mtcars).fit()).to_dataframe().round(3)

Out[4]:

	Term	VIF	VIF_CI_low	VIF_CI_high	SE_factor	Tolerance	Tolerance_CI_low	Tolerance_CI_high
0	wt	2.581	1.778	4.212	1.607	0.387	0.237	0.562
1	hp	3.496	2.308	5.765	1.870	0.286	0.173	0.433
2	C(cyl)	5.106	3.244	8.513	1.503	0.196	0.117	0.308

r2(m)
#> # R2 for Linear Regression
#>        R2: 0.827
#>   adj. R2: 0.815
model_performance(m)

In [5]:

from checkmodel import model_performance, r2

print(r2(m))
model_performance(m)

from checkmodel import model_performance, r2 print(r2(m)) model_performance(m)

# R2 for Linear Regression
       R2: 0.827
  adj. R2: 0.815

Out[5]:

	AIC	AICc	BIC	R2	R2_adjusted	RMSE	Sigma
0	156.652339	158.13382	162.515282	0.826785	0.81484	2.468854	2.593412

m_glm <- glm(vs ~ mpg + wt, family = binomial, data = mtcars)
r2(m_glm)   # Tjur's R2

In [6]:

glm = smf.glm("vs ~ mpg + wt", data=mtcars, family=sm.families.Binomial()).fit()
r2(glm)

glm = smf.glm("vs ~ mpg + wt", data=mtcars, family=sm.families.Binomial()).fit() r2(glm)

Out[6]:

# R2 for Logistic Regression
  Tjur's R2: 0.478

library(lme4)
m_mixed <- lmer(Reaction ~ Days + (Days | Subject), data = sleepstudy)
icc(m_mixed)

In [7]:

from checkmodel import icc

sleep = pd.read_csv(f"{DATA}/sleepstudy.csv")
m_mixed = smf.mixedlm("Reaction ~ Days", sleep, groups="Subject", re_formula="~Days").fit()
icc(m_mixed)

from checkmodel import icc sleep = pd.read_csv(f"{DATA}/sleepstudy.csv") m_mixed = smf.mixedlm("Reaction ~ Days", sleep, groups="Subject", re_formula="~Days").fit() icc(m_mixed)

Out[7]:

# Intraclass Correlation Coefficient
    Adjusted ICC: 0.722
  Unadjusted ICC: 0.521

Where the results differ¶

A handful of defaults differ on purpose: the autocorrelation p-value is exact instead of bootstrapped (and, as in R, check_model() only tests autocorrelation for time-indexed data), a rank-deficient model fails the collinearity check, penalized fits come with a warning, heteroscedasticity of mixed models uses the model's sigma, and a few Bayesian metrics are not ported. See docs/comparison-with-r.md for the full list.