Sesame street

Causal analysis of Sesame Street experiment. See Chapters 18 and 21 in Regression and Other Stories.

Source: Sesame/sesame.Rmd.

The Python version keeps data handling explicit and uses lapylace for Stan-backed generalized linear models, so the statistical model can be read from a formula rather than from handwritten Stan.

from pathlib import Path
import sys
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import lapylace as lp

sys.path.append(str(Path.cwd().parent / 'python'))
from data import ros_path

Data

sesame = pd.read_csv(ros_path('Sesame/data', 'sesame.csv'))
sesame.head()
rownames id site sex age viewcat setting viewenc prebody prelet ... encour _Isite_2 _Isite_3 _Isite_4 _Isite_5 regular watched encouraged y pretest
0 1 1 1 1 66 1 2 1 16 23 ... 1 0 0 0 0 0 0 1 30 23
1 2 2 1 2 67 3 2 1 30 26 ... 1 0 0 0 0 1 1 1 37 26
2 3 3 1 1 56 3 2 2 22 14 ... 0 0 0 0 0 1 1 0 46 14
3 4 4 1 1 49 1 2 2 23 11 ... 0 0 0 0 0 0 0 0 14 11
4 5 5 1 1 69 4 2 2 32 47 ... 0 0 0 0 0 1 1 0 63 47

5 rows × 32 columns

Models

The formulas below are the Python counterparts of the model formulas in the source example. Use lapylace for the Stan-backed Bayesian fit with the same formula interface.

# fit = lp.stan_glm('watched ~ encouraged', data=sesame, family=lp.gaussian(), chains=4, iter_sampling=1000)

# fit = lp.stan_glm('postlet ~ encouraged', data=sesame, family=lp.gaussian(), chains=4, iter_sampling=1000)

# fit = lp.stan_glm('watched ~ encouraged', data=sesame, family=lp.gaussian(), chains=4, iter_sampling=1000)

# fit = lp.stan_glm('postlet ~ watched_hat', data=sesame, family=lp.gaussian(), chains=4, iter_sampling=1000)

# fit = lp.stan_glm('watched ~ encouraged + prelet + as_factor(site) + setting', data=sesame, family=lp.gaussian(), chains=4, iter_sampling=1000)

# fit = lp.stan_glm('postlet ~ watched_hat_3 + prelet + as_factor(site) + setting', data=sesame, family=lp.gaussian(), chains=4, iter_sampling=1000)

Notes

  • Source computation blocks represented: 19.
  • Data paths are expressed through the shared ros_path() helper.
  • Formula-based Bayesian regressions are routed through lapylace.stan_glm().
  • Plotting and simulation work uses NumPy, pandas, matplotlib idioms.