Richard McElreath
January 11, 2023
840

# Statistical Rethinking 2023 - Lecture 04

January 11, 2023

## Transcript

1. Statistical Rethinking
4. Categories & Curves
2023

2. www.3blue1brown.com

3. Drawing Inferences
How to use statistical models to get at
scientific estimands?
Need to incorporate causal thinking into
(1) How we draw the statistical models
(2) How we process the results

4. Drawing Inferences
(1) How we draw the statistical models
Generative model + multiple estimands
=> multiple estimators
(2) How we process the results
Only very simple causal estimates show up
in summary tables

5. Categories & Curves
Linear models can do extra-linear things
Categories (indicator & index variables)
We require these tools to build causal
estimators

6. Categories
How to cope with causes that are
not continuous?
Categories: discrete, unordered
types
Want to stratify by category:
Fit a separate line for each

7. 140 150 160 170 180
30 35 40 45 50 55 60
height (cm)
weight (kg)
women
men

8. Think scientifically first
How are height, weight, and sex causally related?
How are height, weight, and sex statistically related?
W
H
S

9. The causes aren’t in the data
W
H
W
H
140 150 160 170 180
30 35 40 45 50 55 60
height (cm)
weight (kg)

10. The causes aren’t in the data
W
H
W
H
140 150 160 170 180
30 35 40 45 50 55 60
height (cm)
weight (kg)

11. The causes aren’t in the data
S
H
S
H
women
men
140 150 160 170 180
0.00 0.04 0.08
height (cm)
Density

12. The causes aren’t in the data
S
H
S
H
women
men
30 40 50 60
0.00 0.04 0.08
weight (kg)
Density
140 150 160 170 180
0.00 0.04 0.08
height (cm)
Density
women
men
S
W
S
W

13. The causes aren’t in the data
S
H
S
H
women
men
30 40 50 60
0.00 0.04 0.08
weight (kg)
Density
140 150 160 170 180
0.00 0.04 0.08
height (cm)
Density
women
men
S
W
S
W

14. W
H
S
height
influences
weight
sex influences both
height & weight
weight is
influenced
by both
height & sex

15. W
H
S
H = fH
(S)
W = fW
(H, S)

16. height
sex
weight

17. W
H
S
H = fH
(S, U)
W = fW
(H, S, V)
U
V
W
unobserved
causes
ignorable unless shared
S = fS
(W)

18. W
H
S T
V
temperature
U
Figure: Lockley & Eizaguirre 2021

19. W
H
S
H = fH
(S, U)
W = fW
(H, S, V)
U
V
W
S = fS
(W)

20. W
H
S
U
V
W
# S=1 female; S=2 male
sim_HW N H W data.frame(S,H,W)
}
Synthetic people

21. W
H
S
U
V
W
# S=1 female; S=2 male
sim_HW N H W data.frame(S,H,W)
}
Synthetic people
for each S, diﬀerent height-weight line

22. W
H
S
U
V
W
# S=1 female; S=2 male
sim_HW N H W data.frame(S,H,W)
}
S H W
1 1 149.2428 80.07163
2 2 164.8629 94.79647
3 2 158.1800 90.59944
4 1 148.2153 70.68285
5 2 155.5794 91.43798
6 2 154.0733 92.18237

23. Think scientifically first
Diﬀerent causal questions may need
diﬀerent statistical models
Q: Causal eﬀect of H on W?
Q: Causal eﬀect of S on W?
Q: Direct causal eﬀect of S on W?
W
H
S

24. Think scientifically first
Diﬀerent causal questions need
diﬀerent statistical models
Q: Causal eﬀect of H on W?
Q: Causal eﬀect of S on W?
Q: Direct causal eﬀect of S on W?
W
H
S

25. Think scientifically first
Diﬀerent causal questions need
diﬀerent statistical models
Q: Causal eﬀect of H on W?
Q: Causal eﬀect of S on W?
Q: Direct causal eﬀect of S on W?
W
H
S

26. Think scientifically first
Diﬀerent causal questions need
diﬀerent statistical models
Q: Causal eﬀect of H on W?
Q: Causal eﬀect of S on W?
Q: Direct causal eﬀect of S on W?
W
H
S

27. From estimand to estimate
Q: Causal eﬀect of S on W?
W
H
S
W
H
S
Q: Direct causal eﬀect of S on W?
Stratify by S: diﬀerent estimate for each value S can take

28. Drawing Categorical Owls
Several ways to code categorical variables
(1) indicator (0/1) variables
(2) index variables: 1,2,3,4,…
We will use index variables:
Extend to many categories with no change in code
Better for specifying priors
Extend eﬀortlessly to multi-level models

29. Drawing Categorical Owls
Example:
Category Cyan Magenta Yellow Black
Index Value 1 2 3 4
α = [α1
, α2
, α3
, α4
]
Influence of color coded by:

30. Drawing Categorical Owls
α = [α1
, α2
, α3
, α4
]
AAACNHicbVBNS8NAEN34bf2qevSyWBQFKYmIehEEL4IgClaFJoTJdtsu3U3C7kQsoT/Kiz/EiwgeFPHqb3Abc/DrwcKb92aY2RelUhh03SdnZHRsfGJyaroyMzs3v1BdXLo0SaYZb7BEJvo6AsOliHkDBUp+nWoOKpL8KuodDf2rG66NSOIL7Kc8UNCJRVswQCuF1ZN+KOi6b4SiPvJbzE8TrUAONnyVhWLLGh0Fm75fKWq6fkB9kGkXwrxoNywvbhg0RTAIqzW37hagf4lXkhopcRZWH/xWwjLFY2QSjGl6bopBDhoFk3xQ8TPDU2A96PCmpTEoboJyIV2zSou2E21fjLRQv0/koIzpq8h2KsCu+e0Nxf+8Zobt/SAXcZohj9nXonYmKSZ0mCBtCc0Zyr4lwLSwt1LWBQ0Mbc4VG4L3+8t/yeV23dut75zv1A6PyzimyApZJRvEI3vkkByTM9IgjNyRR/JCXp1759l5c96/WkeccmaZ/IDz8Qn0zKu5
yi
⇠ Normal(µi, )
µi = ↵color[i]

31. Using Index Variables
Wi
⇠ Normal(µi, )
µi = ↵ + (Hi
¯
H)
intercept

32. Using Index Variables
AAACQ3icbVBNaxsxENUmbZM4/XCbYy6iJsGhrdktIe0lENKLTyWldRxqLcusPLZFpN1Fmg0xi/9bLvkDueUP9JJDSum1UNlxIR99MPDmvRmkeWmhlaMwvAwWFh89frK0vFJbffrs+Yv6y1eHLi+txI7MdW6PUnCoVYYdUqTxqLAIJtXYTY8/Tf3uCVqn8uwbjQuMDQwzNVASyEtJ/Xs3UXxTOGW4IDyl6nNuDehJU5gyUW+9MTSwJURt1vPNXS5AFyNIqq89FU/4Gy5SpH9ts52od14By9tbSb0RtsIZ+EMSzUmDzXGQ1C9EP5elwYykBud6UVhQXIElJTVOaqJ0WIA8hiH2PM3AoIurWQYTvuGVPh/k1ldGfKbe3qjAODc2qZ80QCN335uK//N6JQ0+xpXKipIwkzcPDUrNKefTQHlfWZSkx56AtMr/lcsRWJDkY6/5EKL7Jz8kh+9b0U5r+8t2Y29/HscyW2evWZNF7APbY212wDpMsjP2g12zn8F5cBX8Cn7fjC4E8501dgfBn78/na9B
Wi
⇠ Normal(µi, )
µi = ↵S[i]
+ S[i]
(Hi
¯
H)
sex of i-th
person
AAACDHicbZDLSsNAFIYn9VbrrerSzWARXEhJSlE3QtGNywr2Ak0oJ9NJO3QyCTMToYQ+gBtfxY0LRdz6AO58GydtFtr6w8DHf85hzvn9mDOlbfvbKqysrq1vFDdLW9s7u3vl/YO2ihJJaItEPJJdHxTlTNCWZprTbiwphD6nHX98k9U7D1QqFol7PYmpF8JQsIAR0MbqlytuCHrkBy7weAT4Cvfm1HfOcqh5psuu2jPhZXByqKBczX75yx1EJAmp0ISDUj3HjrWXgtSMcDotuYmiMZAxDGnPoICQKi+dHTPFJ8YZ4CCS5gmNZ+7viRRCpSahbzqz1dViLTP/q/USHVx6KRNxoqkg84+ChGMd4SwZPGCSEs0nBoBIZnbFZAQSiDb5lUwIzuLJy9CuVZ3zav2uXmlc53EU0RE6RqfIQReogW5RE7UQQY/oGb2iN+vJerHerY95a8HKZw7RH1mfP1Wnmog=
↵ = [↵1, ↵2]
S = 1 indicates female
S = 2 indicates male
Two intercepts, one
for each value in S

33. Using Index Variables
AAACQ3icbVBNaxsxENUmbZM4/XCbYy6iJsGhrdktIe0lENKLTyWldRxqLcusPLZFpN1Fmg0xi/9bLvkDueUP9JJDSum1UNlxIR99MPDmvRmkeWmhlaMwvAwWFh89frK0vFJbffrs+Yv6y1eHLi+txI7MdW6PUnCoVYYdUqTxqLAIJtXYTY8/Tf3uCVqn8uwbjQuMDQwzNVASyEtJ/Xs3UXxTOGW4IDyl6nNuDehJU5gyUW+9MTSwJURt1vPNXS5AFyNIqq89FU/4Gy5SpH9ts52od14By9tbSb0RtsIZ+EMSzUmDzXGQ1C9EP5elwYykBud6UVhQXIElJTVOaqJ0WIA8hiH2PM3AoIurWQYTvuGVPh/k1ldGfKbe3qjAODc2qZ80QCN335uK//N6JQ0+xpXKipIwkzcPDUrNKefTQHlfWZSkx56AtMr/lcsRWJDkY6/5EKL7Jz8kh+9b0U5r+8t2Y29/HscyW2evWZNF7APbY212wDpMsjP2g12zn8F5cBX8Cn7fjC4E8501dgfBn78/na9B
Wi
⇠ Normal(µi, )
µi = ↵S[i]
+ S[i]
(Hi
¯
H)
AAACDHicbZDLSsNAFIYn9VbrrerSzWARXEhJSlE3QtGNywr2Ak0oJ9NJO3QyCTMToYQ+gBtfxY0LRdz6AO58GydtFtr6w8DHf85hzvn9mDOlbfvbKqysrq1vFDdLW9s7u3vl/YO2ihJJaItEPJJdHxTlTNCWZprTbiwphD6nHX98k9U7D1QqFol7PYmpF8JQsIAR0MbqlytuCHrkBy7weAT4Cvfm1HfOcqh5psuu2jPhZXByqKBczX75yx1EJAmp0ISDUj3HjrWXgtSMcDotuYmiMZAxDGnPoICQKi+dHTPFJ8YZ4CCS5gmNZ+7viRRCpSahbzqz1dViLTP/q/USHVx6KRNxoqkg84+ChGMd4SwZPGCSEs0nBoBIZnbFZAQSiDb5lUwIzuLJy9CuVZ3zav2uXmlc53EU0RE6RqfIQReogW5RE7UQQY/oGb2iN+vJerHerY95a8HKZw7RH1mfP1Wnmog=
↵ = [↵1, ↵2]
i = 1
S[1]=2

34. Using Index Variables
AAACQ3icbVBNaxsxENUmbZM4/XCbYy6iJsGhrdktIe0lENKLTyWldRxqLcusPLZFpN1Fmg0xi/9bLvkDueUP9JJDSum1UNlxIR99MPDmvRmkeWmhlaMwvAwWFh89frK0vFJbffrs+Yv6y1eHLi+txI7MdW6PUnCoVYYdUqTxqLAIJtXYTY8/Tf3uCVqn8uwbjQuMDQwzNVASyEtJ/Xs3UXxTOGW4IDyl6nNuDehJU5gyUW+9MTSwJURt1vPNXS5AFyNIqq89FU/4Gy5SpH9ts52od14By9tbSb0RtsIZ+EMSzUmDzXGQ1C9EP5elwYykBud6UVhQXIElJTVOaqJ0WIA8hiH2PM3AoIurWQYTvuGVPh/k1ldGfKbe3qjAODc2qZ80QCN335uK//N6JQ0+xpXKipIwkzcPDUrNKefTQHlfWZSkx56AtMr/lcsRWJDkY6/5EKL7Jz8kh+9b0U5r+8t2Y29/HscyW2evWZNF7APbY212wDpMsjP2g12zn8F5cBX8Cn7fjC4E8501dgfBn78/na9B
Wi
⇠ Normal(µi, )
µi = ↵S[i]
+ S[i]
(Hi
¯
H)
AAACDHicbZDLSsNAFIYn9VbrrerSzWARXEhJSlE3QtGNywr2Ak0oJ9NJO3QyCTMToYQ+gBtfxY0LRdz6AO58GydtFtr6w8DHf85hzvn9mDOlbfvbKqysrq1vFDdLW9s7u3vl/YO2ihJJaItEPJJdHxTlTNCWZprTbiwphD6nHX98k9U7D1QqFol7PYmpF8JQsIAR0MbqlytuCHrkBy7weAT4Cvfm1HfOcqh5psuu2jPhZXByqKBczX75yx1EJAmp0ISDUj3HjrWXgtSMcDotuYmiMZAxDGnPoICQKi+dHTPFJ8YZ4CCS5gmNZ+7viRRCpSahbzqz1dViLTP/q/USHVx6KRNxoqkg84+ChGMd4SwZPGCSEs0nBoBIZnbFZAQSiDb5lUwIzuLJy9CuVZ3zav2uXmlc53EU0RE6RqfIQReogW5RE7UQQY/oGb2iN+vJerHerY95a8HKZw7RH1mfP1Wnmog=
↵ = [↵1, ↵2]
i = 2
S[2]=1

35. Using Index Variables
AAACQ3icbVBNaxsxENUmbZM4/XCbYy6iJsGhrdktIe0lENKLTyWldRxqLcusPLZFpN1Fmg0xi/9bLvkDueUP9JJDSum1UNlxIR99MPDmvRmkeWmhlaMwvAwWFh89frK0vFJbffrs+Yv6y1eHLi+txI7MdW6PUnCoVYYdUqTxqLAIJtXYTY8/Tf3uCVqn8uwbjQuMDQwzNVASyEtJ/Xs3UXxTOGW4IDyl6nNuDehJU5gyUW+9MTSwJURt1vPNXS5AFyNIqq89FU/4Gy5SpH9ts52od14By9tbSb0RtsIZ+EMSzUmDzXGQ1C9EP5elwYykBud6UVhQXIElJTVOaqJ0WIA8hiH2PM3AoIurWQYTvuGVPh/k1ldGfKbe3qjAODc2qZ80QCN335uK//N6JQ0+xpXKipIwkzcPDUrNKefTQHlfWZSkx56AtMr/lcsRWJDkY6/5EKL7Jz8kh+9b0U5r+8t2Y29/HscyW2evWZNF7APbY212wDpMsjP2g12zn8F5cBX8Cn7fjC4E8501dgfBn78/na9B
Wi
⇠ Normal(µi, )
µi = ↵S[i]
+ S[i]
(Hi
¯
H)
Priors
AAACDHicbZDLSsNAFIYn9VbrrerSzWARXEhJSlE3QtGNywr2Ak0oJ9NJO3QyCTMToYQ+gBtfxY0LRdz6AO58GydtFtr6w8DHf85hzvn9mDOlbfvbKqysrq1vFDdLW9s7u3vl/YO2ihJJaItEPJJdHxTlTNCWZprTbiwphD6nHX98k9U7D1QqFol7PYmpF8JQsIAR0MbqlytuCHrkBy7weAT4Cvfm1HfOcqh5psuu2jPhZXByqKBczX75yx1EJAmp0ISDUj3HjrWXgtSMcDotuYmiMZAxDGnPoICQKi+dHTPFJ8YZ4CCS5gmNZ+7viRRCpSahbzqz1dViLTP/q/USHVx6KRNxoqkg84+ChGMd4SwZPGCSEs0nBoBIZnbFZAQSiDb5lUwIzuLJy9CuVZ3zav2uXmlc53EU0RE6RqfIQReogW5RE7UQQY/oGb2iN+vJerHerY95a8HKZw7RH1mfP1Wnmog=
↵ = [↵1, ↵2] AAACC3icbVA9SwNBEN3zM8avqKXNEhUiSLgTUcugjZUoGCPkQpjbbMya3btjd04MR3ob/4qNoCJ24h+w84do7Sax0OiDgcd7M8zMC2IpDLruuzMyOjY+MZmZyk7PzM7N5xYWT02UaMbLLJKRPgvAcClCXkaBkp/FmoMKJK8E7f2eX7nk2ogoPMFOzGsKzkPRFAzQSvVc3gcZt6B+QX0jFPWRX2F6GGkFslvYdjc8d72eW3GLbh/0L/G+yUpp9eP+5XL686iee/MbEUsUD5FJMKbquTHWUtAomOTdrJ8YHgNrwzmvWhqC4qaW9n/p0jWrNGgz0rZCpH3150QKypiOCmynAmyZYa8n/udVE2zu1lIRxgnykA0WNRNJMaK9YGhDaM5QdiwBpoW9lbIWaGBo48vaELzhl/+S082it13cOrZp7JEBMmSZ5EmBeGSHlMgBOSJlwsg1uSUP5NG5ce6cJ+d50DrifM8skV9wXr8A+UieIw==
↵j
⇠ Normal(60, 10)

36. Testing
AAACQ3icbVBNaxsxENUmbZM4/XCbYy6iJsGhrdktIe0lENKLTyWldRxqLcusPLZFpN1Fmg0xi/9bLvkDueUP9JJDSum1UNlxIR99MPDmvRmkeWmhlaMwvAwWFh89frK0vFJbffrs+Yv6y1eHLi+txI7MdW6PUnCoVYYdUqTxqLAIJtXYTY8/Tf3uCVqn8uwbjQuMDQwzNVASyEtJ/Xs3UXxTOGW4IDyl6nNuDehJU5gyUW+9MTSwJURt1vPNXS5AFyNIqq89FU/4Gy5SpH9ts52od14By9tbSb0RtsIZ+EMSzUmDzXGQ1C9EP5elwYykBud6UVhQXIElJTVOaqJ0WIA8hiH2PM3AoIurWQYTvuGVPh/k1ldGfKbe3qjAODc2qZ80QCN335uK//N6JQ0+xpXKipIwkzcPDUrNKefTQHlfWZSkx56AtMr/lcsRWJDkY6/5EKL7Jz8kh+9b0U5r+8t2Y29/HscyW2evWZNF7APbY212wDpMsjP2g12zn8F5cBX8Cn7fjC4E8501dgfBn78/na9B
Wi
⇠ Normal(µi, )
µi = ↵S[i]
+ S[i]
(Hi
¯
H)
AAACCnicbVDLSgMxFM3UV62vUZduokWpIGVGirosunFZwWkLnaFk0rQNTWaG5I5Yhq7d+CtuXCji1i9w59+YPhZavRByOOdc7r0nTATX4DhfVm5hcWl5Jb9aWFvf2Nyyt3fqOk4VZR6NRayaIdFM8Ih5wEGwZqIYkaFgjXBwNdYbd0xpHke3MExYIEkv4l1OCRiqbe/7mvckwUfml9gHdg+ZZwyxkqOSc+I6x2276JSdSeG/wJ2BIppVrW1/+p2YppJFQAXRuuU6CQQZUcCpYKOCn2qWEDogPdYyMCKS6SCbnDLCh4bpYDPevAjwhP3ZkRGp9VCGxikJ9PW8Nib/01opdC+CjEdJCiyi00HdVGCI8TgX3OGKURBDAwhV3OyKaZ8oQsGkVzAhuPMn/wX107J7Vq7cVIrVy1kcebSHDlAJuegcVdE1qiEPUfSAntALerUerWfrzXqfWnPWrGcX/Srr4xt9+5l9
⇠ Uniform(0, 10)
AAACC3icbVA9SwNBEN3zM8avqKXNEhUiSLgTUcugjZUoGCPkQpjbbMya3btjd04MR3ob/4qNoCJ24h+w84do7Sax0OiDgcd7M8zMC2IpDLruuzMyOjY+MZmZyk7PzM7N5xYWT02UaMbLLJKRPgvAcClCXkaBkp/FmoMKJK8E7f2eX7nk2ogoPMFOzGsKzkPRFAzQSvVc3gcZt6B+QX0jFPWRX2F6GGkFslvYdjc8d72eW3GLbh/0L/G+yUpp9eP+5XL686iee/MbEUsUD5FJMKbquTHWUtAomOTdrJ8YHgNrwzmvWhqC4qaW9n/p0jWrNGgz0rZCpH3150QKypiOCmynAmyZYa8n/udVE2zu1lIRxgnykA0WNRNJMaK9YGhDaM5QdiwBpoW9lbIWaGBo48vaELzhl/+S082it13cOrZp7JEBMmSZ5EmBeGSHlMgBOSJlwsg1uSUP5NG5ce6cJ+d50DrifM8skV9wXr8A+UieIw==
↵j
⇠ Normal(60, 10)
# female sample
S simF # male sample
S simM # effect of sex (male-female)
mean(simM\$W - simF\$W)
What is the total causal eﬀect of sex? Causal
eﬀect of sex is diﬀerence made intervening:
[1] 21.13852

37. # observe sample
S dat # estimate posterior
m_SW alist(
W ~ dnorm(mu,sigma),
mu a[S] ~ dnorm(60,10),
sigma ~ dunif(0,10)
), data=dat )
precis(m_SW,depth=2)
Now run the estimating model and synthetic sample
a[1] 75.47 0.93 73.98 76.95
a[2] 97.39 0.84 96.05 98.74
sigma 6.26 0.44 5.55 6.96
AAACQ3icbVBNaxsxENUmbZM4/XCbYy6iJsGhrdktIe0lENKLTyWldRxqLcusPLZFpN1Fmg0xi/9bLvkDueUP9JJDSum1UNlxIR99MPDmvRmkeWmhlaMwvAwWFh89frK0vFJbffrs+Yv6y1eHLi+txI7MdW6PUnCoVYYdUqTxqLAIJtXYTY8/Tf3uCVqn8uwbjQuMDQwzNVASyEtJ/Xs3UXxTOGW4IDyl6nNuDehJU5gyUW+9MTSwJURt1vPNXS5AFyNIqq89FU/4Gy5SpH9ts52od14By9tbSb0RtsIZ+EMSzUmDzXGQ1C9EP5elwYykBud6UVhQXIElJTVOaqJ0WIA8hiH2PM3AoIurWQYTvuGVPh/k1ldGfKbe3qjAODc2qZ80QCN335uK//N6JQ0+xpXKipIwkzcPDUrNKefTQHlfWZSkx56AtMr/lcsRWJDkY6/5EKL7Jz8kh+9b0U5r+8t2Y29/HscyW2evWZNF7APbY212wDpMsjP2g12zn8F5cBX8Cn7fjC4E8501dgfBn78/na9B
Wi
⇠ Normal(µi, )
µi = ↵S[i]
+ S[i]
(Hi
¯
H)
AAACCnicbVDLSgMxFM3UV62vUZduokWpIGVGirosunFZwWkLnaFk0rQNTWaG5I5Yhq7d+CtuXCji1i9w59+YPhZavRByOOdc7r0nTATX4DhfVm5hcWl5Jb9aWFvf2Nyyt3fqOk4VZR6NRayaIdFM8Ih5wEGwZqIYkaFgjXBwNdYbd0xpHke3MExYIEkv4l1OCRiqbe/7mvckwUfml9gHdg+ZZwyxkqOSc+I6x2276JSdSeG/wJ2BIppVrW1/+p2YppJFQAXRuuU6CQQZUcCpYKOCn2qWEDogPdYyMCKS6SCbnDLCh4bpYDPevAjwhP3ZkRGp9VCGxikJ9PW8Nib/01opdC+CjEdJCiyi00HdVGCI8TgX3OGKURBDAwhV3OyKaZ8oQsGkVzAhuPMn/wX107J7Vq7cVIrVy1kcebSHDlAJuegcVdE1qiEPUfSAntALerUerWfrzXqfWnPWrGcX/Srr4xt9+5l9
⇠ Uniform(0, 10)
AAACC3icbVA9SwNBEN3zM8avqKXNEhUiSLgTUcugjZUoGCPkQpjbbMya3btjd04MR3ob/4qNoCJ24h+w84do7Sax0OiDgcd7M8zMC2IpDLruuzMyOjY+MZmZyk7PzM7N5xYWT02UaMbLLJKRPgvAcClCXkaBkp/FmoMKJK8E7f2eX7nk2ogoPMFOzGsKzkPRFAzQSvVc3gcZt6B+QX0jFPWRX2F6GGkFslvYdjc8d72eW3GLbh/0L/G+yUpp9eP+5XL686iee/MbEUsUD5FJMKbquTHWUtAomOTdrJ8YHgNrwzmvWhqC4qaW9n/p0jWrNGgz0rZCpH3150QKypiOCmynAmyZYa8n/udVE2zu1lIRxgnykA0WNRNJMaK9YGhDaM5QdiwBpoW9lbIWaGBo48vaELzhl/+S082it13cOrZp7JEBMmSZ5EmBeGSHlMgBOSJlwsg1uSUP5NG5ce6cJ+d50DrifM8skV9wXr8A+UieIw==
↵j
⇠ Normal(60, 10)
97 – 75 = 22

38. Analyze the real sample
AAACQ3icbVBNaxsxENUmbZM4/XCbYy6iJsGhrdktIe0lENKLTyWldRxqLcusPLZFpN1Fmg0xi/9bLvkDueUP9JJDSum1UNlxIR99MPDmvRmkeWmhlaMwvAwWFh89frK0vFJbffrs+Yv6y1eHLi+txI7MdW6PUnCoVYYdUqTxqLAIJtXYTY8/Tf3uCVqn8uwbjQuMDQwzNVASyEtJ/Xs3UXxTOGW4IDyl6nNuDehJU5gyUW+9MTSwJURt1vPNXS5AFyNIqq89FU/4Gy5SpH9ts52od14By9tbSb0RtsIZ+EMSzUmDzXGQ1C9EP5elwYykBud6UVhQXIElJTVOaqJ0WIA8hiH2PM3AoIurWQYTvuGVPh/k1ldGfKbe3qjAODc2qZ80QCN335uK//N6JQ0+xpXKipIwkzcPDUrNKefTQHlfWZSkx56AtMr/lcsRWJDkY6/5EKL7Jz8kh+9b0U5r+8t2Y29/HscyW2evWZNF7APbY212wDpMsjP2g12zn8F5cBX8Cn7fjC4E8501dgfBn78/na9B
Wi
⇠ Normal(µi, )
µi = ↵S[i]
+ S[i]
(Hi
¯
H)
data(Howell1)
d d =18 , ]
dat W = d\$weight,
S = d\$male + 1 ) # S=1 female, S=2 male
m_SW alist(
W ~ dnorm(mu,sigma),
mu a[S] ~ dnorm(60,10),
sigma ~ dunif(0,10)
), data=dat ) AAACCnicbVDLSgMxFM3UV62vUZduokWpIGVGirosunFZwWkLnaFk0rQNTWaG5I5Yhq7d+CtuXCji1i9w59+YPhZavRByOOdc7r0nTATX4DhfVm5hcWl5Jb9aWFvf2Nyyt3fqOk4VZR6NRayaIdFM8Ih5wEGwZqIYkaFgjXBwNdYbd0xpHke3MExYIEkv4l1OCRiqbe/7mvckwUfml9gHdg+ZZwyxkqOSc+I6x2276JSdSeG/wJ2BIppVrW1/+p2YppJFQAXRuuU6CQQZUcCpYKOCn2qWEDogPdYyMCKS6SCbnDLCh4bpYDPevAjwhP3ZkRGp9VCGxikJ9PW8Nib/01opdC+CjEdJCiyi00HdVGCI8TgX3OGKURBDAwhV3OyKaZ8oQsGkVzAhuPMn/wX107J7Vq7cVIrVy1kcebSHDlAJuegcVdE1qiEPUfSAntALerUerWfrzXqfWnPWrGcX/Srr4xt9+5l9
⇠ Uniform(0, 10)
AAACC3icbVA9SwNBEN3zM8avqKXNEhUiSLgTUcugjZUoGCPkQpjbbMya3btjd04MR3ob/4qNoCJ24h+w84do7Sax0OiDgcd7M8zMC2IpDLruuzMyOjY+MZmZyk7PzM7N5xYWT02UaMbLLJKRPgvAcClCXkaBkp/FmoMKJK8E7f2eX7nk2ogoPMFOzGsKzkPRFAzQSvVc3gcZt6B+QX0jFPWRX2F6GGkFslvYdjc8d72eW3GLbh/0L/G+yUpp9eP+5XL686iee/MbEUsUD5FJMKbquTHWUtAomOTdrJ8YHgNrwzmvWhqC4qaW9n/p0jWrNGgz0rZCpH3150QKypiOCmynAmyZYa8n/udVE2zu1lIRxgnykA0WNRNJMaK9YGhDaM5QdiwBpoW9lbIWaGBo48vaELzhl/+S082it13cOrZp7JEBMmSZ5EmBeGSHlMgBOSJlwsg1uSUP5NG5ce6cJ+d50DrifM8skV9wXr8A+UieIw==
↵j
⇠ Normal(60, 10)

39. Posterior means & predictions
# posterior mean W
post dens( post\$a[,1] , xlim=c(39,50) , lwd=3 ,
col=2 , xlab="posterior mean weight (kg)" )
dens( post\$a[,2] , lwd=3 , col=4 , add=TRUE )
# posterior W distributions
W1 W2 dens( W1 , xlim=c(20,70) , ylim=c(0,0.085) ,
lwd=3 , col=2 )
dens( W2 , lwd=3 , col=4 , add=TRUE )
40 42 44 46 48 50
0.0 0.4 0.8
posterior mean weight (kg)
Density
women men

40. Posterior means & predictions
# posterior mean W
post dens( post\$a[,1] , xlim=c(39,50) , lwd=3 ,
col=2 , xlab="posterior mean weight (kg)" )
dens( post\$a[,2] , lwd=3 , col=4 , add=TRUE )
# posterior W distributions
W1 W2 dens( W1 , xlim=c(20,70) , ylim=c(0,0.085) ,
lwd=3 , col=2 )
dens( W2 , lwd=3 , col=4 , add=TRUE )
40 42 44 46 48 50
0.0 0.4 0.8
posterior mean weight (kg)
Density
20 30 40 50 60 70
0.00 0.04 0.08
posterior predicted weight (kg)
Density
women men
men
women

41. Always Be Contrasting
Need to compute contrast, the
diﬀerence between the categories
It is never legitimate to compare
overlap in distributions
Must compute contrast
distribution
20 30 40 50 60 70
0.00 0.04 0.08
posterior predicted weight (kg)
Density

42. -1.5 -0.5 0.0 0.5 1.0 1.5
-0.5 0.5 1.0 1.5 2.0 2.5
parameter 1
parameter 2
It is never legitimate to compare overlap in parameters or lines

43. -3 -2 -1 0 1 2 3
0.0 0.5 1.0 1.5 2.0
value
Density
-1.5 -0.5 0.0 0.5 1.0 1.5
-0.5 0.5 1.0 1.5 2.0 2.5
parameter 1
parameter 2
parameter 1
parameter 2
diﬀerence
It is never legitimate to compare overlap in parameters or lines
This means no comparing confidence intervals or p-values either!

44. Causal contrast
# causal contrast (in means)
mu_contrast dens( mu_contrast , xlim=c(3,10) , lwd=3 ,
col=1 , xlab="posterior mean weight contrast
(kg)" )
40 42 44 46 48 50
0.0 0.4 0.8
posterior mean weight (kg)
Density
women men
3 4 5 6 7 8 9 10
0.0 0.3 0.6
posterior mean weight contrast (kg)
Density

45. Weight contrast
20 30 40 50 60 70
0.00 0.04 0.08
posterior predicted weight (kg)
Density
men
women

46. Weight contrast
# posterior W distributions
W1 W2 # contrast
W_contrast dens( W_contrast , xlim=c(-25,35) , lwd=3 ,
col=1 , xlab="posterior weight contrast (kg)" )
# proportion above zero
sum( W_contrast > 0 ) / 1000
# proportion below zero
sum( W_contrast < 0 ) / 1000
20 30 40 50 60 70
0.00 0.04 0.08
posterior predicted weight (kg)
Density
men
women
-20 -10 0 10 20 30
0.00 0.02 0.04
posterior weight contrast (kg)
Density
18% 82%

47. From estimand to estimate
Q: Causal eﬀect of S on W?
W
H
S
W
H
S
Q: Direct causal eﬀect of S on W?
3 4 5 6 7 8 9 10
0.0 0.3 0.6
posterior mean weight contrast (kg)
Density
-20 -10 0 10 20 30
0.00 0.02 0.04
posterior weight contrast (kg)
Density

48. From estimand to estimate
W
H
S
Q: Direct causal eﬀect of S on W?
Now we need to “block”
association through H
This means stratify by H

49. S dat 140 145 150 155 160 165 170
70 80 90 100
H
W
direct
indirect
women
men

50. 140 145 150 155 160 165 170
70 80 90 100
H
W
S dat direct
indirect
Indirect: H influences W
diﬀerently (slopes)
Direct: W diﬀers net any
slope diﬀerences

51. Index Variables & Lines
Wi
⇠ Normal(µi, )
µi = ↵ + (Hi
¯
H)
intercept slope

52. Centering variables
Wi
⇠ Normal(µi, )
µi = ↵ + (Hi
¯
H)
“centered” height
Centering H makes it so that
alpha is the average weight of
a person with average height
average height
Centering makes it easier to
define scientific priors for alpha

53. 140 145 150 155 160 165 170
70 80 90 100
H
W
Centering H makes it so
that alpha is the average
weight of a person with
average height
average W
average H
Wi
⇠ Normal(µi, )
µi = ↵ + (Hi
¯
H)

54. 140 145 150 155 160 165 170
70 80 90 100
H
W
Centering H makes it so
that alpha is the average
weight of a person with
average height
average W
average H
Wi
⇠ Normal(µi, )
µi = ↵ + (Hi
¯
H)

55. Index Variables & Lines
AAACQ3icbVBNaxsxENUmbZM4/XCbYy6iJsGhrdktIe0lENKLTyWldRxqLcusPLZFpN1Fmg0xi/9bLvkDueUP9JJDSum1UNlxIR99MPDmvRmkeWmhlaMwvAwWFh89frK0vFJbffrs+Yv6y1eHLi+txI7MdW6PUnCoVYYdUqTxqLAIJtXYTY8/Tf3uCVqn8uwbjQuMDQwzNVASyEtJ/Xs3UXxTOGW4IDyl6nNuDehJU5gyUW+9MTSwJURt1vPNXS5AFyNIqq89FU/4Gy5SpH9ts52od14By9tbSb0RtsIZ+EMSzUmDzXGQ1C9EP5elwYykBud6UVhQXIElJTVOaqJ0WIA8hiH2PM3AoIurWQYTvuGVPh/k1ldGfKbe3qjAODc2qZ80QCN335uK//N6JQ0+xpXKipIwkzcPDUrNKefTQHlfWZSkx56AtMr/lcsRWJDkY6/5EKL7Jz8kh+9b0U5r+8t2Y29/HscyW2evWZNF7APbY212wDpMsjP2g12zn8F5cBX8Cn7fjC4E8501dgfBn78/na9B
Wi
⇠ Normal(µi, )
µi = ↵S[i]
+ S[i]
(Hi
¯
H)
sex of i-th
person
= [ 1, 2]
AAACDHicbZDLSsNAFIYn9VbrrerSzWARXEhJSlE3QtGNywr2Ak0oJ9NJO3QyCTMToYQ+gBtfxY0LRdz6AO58GydtFtr6w8DHf85hzvn9mDOlbfvbKqysrq1vFDdLW9s7u3vl/YO2ihJJaItEPJJdHxTlTNCWZprTbiwphD6nHX98k9U7D1QqFol7PYmpF8JQsIAR0MbqlytuCHrkBy7weAT4Cvfm1HfOcqh5psuu2jPhZXByqKBczX75yx1EJAmp0ISDUj3HjrWXgtSMcDotuYmiMZAxDGnPoICQKi+dHTPFJ8YZ4CCS5gmNZ+7viRRCpSahbzqz1dViLTP/q/USHVx6KRNxoqkg84+ChGMd4SwZPGCSEs0nBoBIZnbFZAQSiDb5lUwIzuLJy9CuVZ3zav2uXmlc53EU0RE6RqfIQReogW5RE7UQQY/oGb2iN+vJerHerY95a8HKZw7RH1mfP1Wnmog=
↵ = [↵1, ↵2]
S = 1 indicates female
S = 2 indicates male
Two intercepts and two slopes, one for each value in S

56. Index Variables & Lines
AAACQ3icbVBNaxsxENUmbZM4/XCbYy6iJsGhrdktIe0lENKLTyWldRxqLcusPLZFpN1Fmg0xi/9bLvkDueUP9JJDSum1UNlxIR99MPDmvRmkeWmhlaMwvAwWFh89frK0vFJbffrs+Yv6y1eHLi+txI7MdW6PUnCoVYYdUqTxqLAIJtXYTY8/Tf3uCVqn8uwbjQuMDQwzNVASyEtJ/Xs3UXxTOGW4IDyl6nNuDehJU5gyUW+9MTSwJURt1vPNXS5AFyNIqq89FU/4Gy5SpH9ts52od14By9tbSb0RtsIZ+EMSzUmDzXGQ1C9EP5elwYykBud6UVhQXIElJTVOaqJ0WIA8hiH2PM3AoIurWQYTvuGVPh/k1ldGfKbe3qjAODc2qZ80QCN335uK//N6JQ0+xpXKipIwkzcPDUrNKefTQHlfWZSkx56AtMr/lcsRWJDkY6/5EKL7Jz8kh+9b0U5r+8t2Y29/HscyW2evWZNF7APbY212wDpMsjP2g12zn8F5cBX8Cn7fjC4E8501dgfBn78/na9B
Wi
⇠ Normal(µi, )
µi = ↵S[i]
+ S[i]
(Hi
¯
H)
= [ 1, 2]
AAACDHicbZDLSsNAFIYn9VbrrerSzWARXEhJSlE3QtGNywr2Ak0oJ9NJO3QyCTMToYQ+gBtfxY0LRdz6AO58GydtFtr6w8DHf85hzvn9mDOlbfvbKqysrq1vFDdLW9s7u3vl/YO2ihJJaItEPJJdHxTlTNCWZprTbiwphD6nHX98k9U7D1QqFol7PYmpF8JQsIAR0MbqlytuCHrkBy7weAT4Cvfm1HfOcqh5psuu2jPhZXByqKBczX75yx1EJAmp0ISDUj3HjrWXgtSMcDotuYmiMZAxDGnPoICQKi+dHTPFJ8YZ4CCS5gmNZ+7viRRCpSahbzqz1dViLTP/q/USHVx6KRNxoqkg84+ChGMd4SwZPGCSEs0nBoBIZnbFZAQSiDb5lUwIzuLJy9CuVZ3zav2uXmlc53EU0RE6RqfIQReogW5RE7UQQY/oGb2iN+vJerHerY95a8HKZw7RH1mfP1Wnmog=
↵ = [↵1, ↵2]
i = 1
S[1]=2

57. Index Variables & Lines
AAACQ3icbVBNaxsxENUmbZM4/XCbYy6iJsGhrdktIe0lENKLTyWldRxqLcusPLZFpN1Fmg0xi/9bLvkDueUP9JJDSum1UNlxIR99MPDmvRmkeWmhlaMwvAwWFh89frK0vFJbffrs+Yv6y1eHLi+txI7MdW6PUnCoVYYdUqTxqLAIJtXYTY8/Tf3uCVqn8uwbjQuMDQwzNVASyEtJ/Xs3UXxTOGW4IDyl6nNuDehJU5gyUW+9MTSwJURt1vPNXS5AFyNIqq89FU/4Gy5SpH9ts52od14By9tbSb0RtsIZ+EMSzUmDzXGQ1C9EP5elwYykBud6UVhQXIElJTVOaqJ0WIA8hiH2PM3AoIurWQYTvuGVPh/k1ldGfKbe3qjAODc2qZ80QCN335uK//N6JQ0+xpXKipIwkzcPDUrNKefTQHlfWZSkx56AtMr/lcsRWJDkY6/5EKL7Jz8kh+9b0U5r+8t2Y29/HscyW2evWZNF7APbY212wDpMsjP2g12zn8F5cBX8Cn7fjC4E8501dgfBn78/na9B
Wi
⇠ Normal(µi, )
µi = ↵S[i]
+ S[i]
(Hi
¯
H)
= [ 1, 2]
AAACDHicbZDLSsNAFIYn9VbrrerSzWARXEhJSlE3QtGNywr2Ak0oJ9NJO3QyCTMToYQ+gBtfxY0LRdz6AO58GydtFtr6w8DHf85hzvn9mDOlbfvbKqysrq1vFDdLW9s7u3vl/YO2ihJJaItEPJJdHxTlTNCWZprTbiwphD6nHX98k9U7D1QqFol7PYmpF8JQsIAR0MbqlytuCHrkBy7weAT4Cvfm1HfOcqh5psuu2jPhZXByqKBczX75yx1EJAmp0ISDUj3HjrWXgtSMcDotuYmiMZAxDGnPoICQKi+dHTPFJ8YZ4CCS5gmNZ+7viRRCpSahbzqz1dViLTP/q/USHVx6KRNxoqkg84+ChGMd4SwZPGCSEs0nBoBIZnbFZAQSiDb5lUwIzuLJy9CuVZ3zav2uXmlc53EU0RE6RqfIQReogW5RE7UQQY/oGb2iN+vJerHerY95a8HKZw7RH1mfP1Wnmog=
↵ = [↵1, ↵2]
i = 2
S[2]=1

58. Analyze the sample
data(Howell1)
d d =18 , ]
dat W = d\$weight,
H = d\$height,
Hbar = mean(d\$height),
S = d\$male + 1 ) # S=1 female, S=2 male
m_SHW alist(
W ~ dnorm(mu,sigma),
mu a[S] ~ dnorm(60,10),
b[S] ~ dunif(0,1),
sigma ~ dunif(0,10)
), data=dat )
140 150 160 170 180
30 35 40 45 50 55 60
height (cm)
weight (kg)

59. Contrasts at each height
(1) Compute posterior predictive
for women
(2) Compute posterior predictive
for men
(3) Subtract (2) from (1)
(4) Plot distribution at each height
(on right)
140 150 160 170 180
30 35 40 45 50 55 60
height (cm)
weight (kg)

60. Contrasts at each height
130 140 150 160 170 180 190
-6 -4 -2 0 2 4 6 8
height (cm)
weight contrast (F–M)
lines( xseq , apply(muF,2,mean) , lwd=3 , col=2 )
lines( xseq , apply(muM,2,mean) , lwd=3 , col=4 )
mu_contrast plot( NULL , xlim=range(xseq) , ylim=c(-6,8) , xlab="height (cm)"
, ylab="weight contrast (F–M)" )
for ( p in c(0.5,0.6,0.7,0.8,0.9,0.99) )
abline(h=0,lty=2)
women heavier
men heavier
Nearly all of the causal eﬀect of S acts through H

61. From estimand to estimate
Q: Causal eﬀect of S on W?
W
H
S
W
H
S
Q: Direct causal eﬀect of S on W?
3 4 5 6 7 8 9 10
0.0 0.3 0.6
posterior mean weight contrast (kg)
Density
-20 -10 0 10 20 30
0.00 0.02 0.04
posterior weight contrast (kg)
Density
130 140 150 160 170 180 190
-6 -4 -2 0 2 4 6 8
height (cm)
weight contrast (F–M)

62. Categorical variables
Easy to use with index coding
Use samples to compute relevant
contrasts
Always summarize (mean, interval)
as the last step
Want mean diﬀerence and not
diﬀerence of means
140 150 160 170 180
30 35 40 45 50 55 60
height (cm)
weight (kg)

63. PAUSE

64. Curves from lines
H –> W obviously not linear
Linear models can easily fit
curves
But this is not mechanistic
library(rethinking)
data(Howell1)
60 80 100 140 180
10 20 30 40 50 60
height (cm)
weight (kg)

65. Curves from lines
Linear models can easily fit
curves
Two popular strategies
(1) polynomials — awful
(2) splines and generalized
library(rethinking)
data(Howell1)
60 80 100 140 180
10 20 30 40 50 60
height (cm)
weight (kg)

66. Polynomial linear models
Strategy: polynomial functions
Problems: strange symmetries,
explosive uncertainty at edges
No local smoothing; only global
DO NOT USE
μi
= α + β1
xi
+ β2
x2
i

67. 50 100 150 200
0 20 40 60
height (cm)
weight (kg)
μi
= α + β1
Hi
+ β2
H2
i

68. 50 100 150 200
0 20 40 60
height (cm)
weight (kg)
50 100 150 200
0 20 40 60
height (cm)
weight (kg)
μi
= α + β1
Hi
+ β2
H2
i
μi
= α + β1
Hi
+ β2
H2
i
= + b3
H3 + b4
H4

69. (&0.&53*\$ 1&01-&
h
r
V = πr2h
'ĶĴłĿĲ ƉƎƉ ćF i
IVNBO XFJHIU BT B G
SVWJBO .BO XFSF B D
IJT XFJHIU CZ DBMDV
GVODUJPO PG IJT IFJHI
Q PG IFJHIU ćJT NFBOT S = QI 4VCTUJUVUJOH UIJT JOUP UIF GPSNVMB
7 = π(QI)I = πQI
'JOBMMZ XFJHIU JT TPNF QSPQPSUJPO PG WPMVNFIPX NBOZ LJMPHSBN
UJNFUFS 4P XF OFFE B QBSBNFUFS L UP FYQSFTT UIJT USBOTMBUJPO CFUX
8 = L7 = LπQI
"OE UIJT JT PVS GPSNVMB GPS FYQFDUFE XFJHIU HJWFO BO JOEJWJEVBMT
WJPVTMZ BO PSEJOBSZ HFOFSBMJ[FE MJOFBS NPEFM #VU UIBUT PLBZ *U
Thinking vs Fitting
Linear models can fit anything
(Geocentrism)
Better to think
60 80 100 140 180
10 20 30 40 50 60
height (cm)
weight (kg)
log Wi
∼ Normal(μi
, σ)
μi
= α + β(Hi
− ¯
H)
Will revisit in lecture 19

70. Splines
Basis-Splines: Wiggly function
built from many local
functions
Basis function: A local
function
Local functions make splines
better than polynomials, but
equally geocentric
Richard McElreath
McElreath
with Examples in R and Stan
SECOND EDITION
ical Rethinking
d Gaussian process models
h (DAG) approach to causal
w edition also contains new
ered categorical predictors,
ampling, instrumental vari-
rely new chapter that goes
n-specific scientific models
mples
umptions are reflected in
in optional sections
analyze causal graphs
bsite and on GitHub
nd is a Director at the Max
Germany. He has published
analysis of social behavior,
Models of Social Evolution.

71. Going Local — B-Splines
B-Splines are just linear models, but with some weird
synthetic variables:
Weights w are like slopes
Basis functions B are synthetic variables
B values turn on weights in diﬀerent regions
Detailed example starting on page 114
CZ HFOFSBUJOH OFX QSFEJDUPS WBSJBCMFT BOE VTJOH UIPTF JO UIF MJOFBS NPEFM µJ
6O
PNJBM SFHSFTTJPO #TQMJOFT EP OPU EJSFDUMZ USBOTGPSN UIF QSFEJDUPS CZ TRVBSJOH PS
*OTUFBE UIFZ JOWFOU B TFSJFT PG FOUJSFMZ OFX TZOUIFUJD QSFEJDUPS WBSJBCMFT &BDI
BSJBCMFT TFSWFT UP HSBEVBMMZ UVSO B TQFDJĕD QBSBNFUFS PO BOE PČ XJUIJO B TQFDJĕD
IF QSFEJDUPS WBSJBCMFT &BDI PG UIFTF WBSJBCMFT JT DBMMFE B įĮŀĶŀ ĳłĻİŁĶļĻ ćF
EFM FOET VQ MPPLJOH WFSZ GBNJMJBS
µJ = α + X
#J, + X
#J, + X
#J, + ...
JT UIF OUI CBTJT GVODUJPOT WBMVF PO SPX J BOE UIF X QBSBNFUFST BSF DPSSFTQPOE
T GPS FBDI ćF QBSBNFUFST BDU MJLF TMPQFT BEKVTUJOH UIF JOĘVFODF PG FBDI CBTJT
O UIF NFBO µJ
4P SFBMMZ UIJT JT KVTU BOPUIFS MJOFBS SFHSFTTJPO CVU XJUI TPNF GBODZ
QSFEJDUPS WBSJBCMFT ćFTF TZOUIFUJD WBSJBCMFT EP TPNF SFBMMZ FMFHBOU EFTDSJQUJWF
XPSL GPS VT
EP XF DPOTUSVDU UIFTF CBTJT WBSJBCMFT # * EJTQMBZ UIF TJNQMFTU DBTF JO 'ĶĴłĿĲ ƌƉƊ
BQQSPYJNBUF UIF UFNQFSBUVSF EBUB XJUI GPVS EJČFSFOU MJOFBS BQQSPYJNBUJPOT 'JSTU

72. Going Local — B-Splines
B-Splines are just linear models, but with some weird
synthetic variables:
Weights w are like slopes
Basis functions B are synthetic variables
B values turn on weights in diﬀerent regions
Detailed example starting on page 114
CZ HFOFSBUJOH OFX QSFEJDUPS WBSJBCMFT BOE VTJOH UIPTF JO UIF MJOFBS NPEFM µJ
6O
PNJBM SFHSFTTJPO #TQMJOFT EP OPU EJSFDUMZ USBOTGPSN UIF QSFEJDUPS CZ TRVBSJOH PS
*OTUFBE UIFZ JOWFOU B TFSJFT PG FOUJSFMZ OFX TZOUIFUJD QSFEJDUPS WBSJBCMFT &BDI
BSJBCMFT TFSWFT UP HSBEVBMMZ UVSO B TQFDJĕD QBSBNFUFS PO BOE PČ XJUIJO B TQFDJĕD
IF QSFEJDUPS WBSJBCMFT &BDI PG UIFTF WBSJBCMFT JT DBMMFE B įĮŀĶŀ ĳłĻİŁĶļĻ ćF
EFM FOET VQ MPPLJOH WFSZ GBNJMJBS
µJ = α + X
#J, + X
#J, + X
#J, + ...
JT UIF OUI CBTJT GVODUJPOT WBMVF PO SPX J BOE UIF X QBSBNFUFST BSF DPSSFTQPOE
T GPS FBDI ćF QBSBNFUFST BDU MJLF TMPQFT BEKVTUJOH UIF JOĘVFODF PG FBDI CBTJT
O UIF NFBO µJ
4P SFBMMZ UIJT JT KVTU BOPUIFS MJOFBS SFHSFTTJPO CVU XJUI TPNF GBODZ
QSFEJDUPS WBSJBCMFT ćFTF TZOUIFUJD WBSJBCMFT EP TPNF SFBMMZ FMFHBOU EFTDSJQUJWF
XPSL GPS VT
EP XF DPOTUSVDU UIFTF CBTJT WBSJBCMFT # * EJTQMBZ UIF TJNQMFTU DBTF JO 'ĶĴłĿĲ ƌƉƊ
BQQSPYJNBUF UIF UFNQFSBUVSF EBUB XJUI GPVS EJČFSFOU MJOFBS BQQSPYJNBUJPOT 'JSTU

73. Going Local — B-Splines
B-Splines are just linear models, but with some weird
synthetic variables:
Weights w are like slopes
Basis functions B are synthetic variables
B values turn on weights in diﬀerent regions
Detailed example starting on page 114
CZ HFOFSBUJOH OFX QSFEJDUPS WBSJBCMFT BOE VTJOH UIPTF JO UIF MJOFBS NPEFM µJ
6O
PNJBM SFHSFTTJPO #TQMJOFT EP OPU EJSFDUMZ USBOTGPSN UIF QSFEJDUPS CZ TRVBSJOH PS
*OTUFBE UIFZ JOWFOU B TFSJFT PG FOUJSFMZ OFX TZOUIFUJD QSFEJDUPS WBSJBCMFT &BDI
BSJBCMFT TFSWFT UP HSBEVBMMZ UVSO B TQFDJĕD QBSBNFUFS PO BOE PČ XJUIJO B TQFDJĕD
IF QSFEJDUPS WBSJBCMFT &BDI PG UIFTF WBSJBCMFT JT DBMMFE B įĮŀĶŀ ĳłĻİŁĶļĻ ćF
EFM FOET VQ MPPLJOH WFSZ GBNJMJBS
µJ = α + X
#J, + X
#J, + X
#J, + ...
JT UIF OUI CBTJT GVODUJPOT WBMVF PO SPX J BOE UIF X QBSBNFUFST BSF DPSSFTQPOE
T GPS FBDI ćF QBSBNFUFST BDU MJLF TMPQFT BEKVTUJOH UIF JOĘVFODF PG FBDI CBTJT
O UIF NFBO µJ
4P SFBMMZ UIJT JT KVTU BOPUIFS MJOFBS SFHSFTTJPO CVU XJUI TPNF GBODZ
QSFEJDUPS WBSJBCMFT ćFTF TZOUIFUJD WBSJBCMFT EP TPNF SFBMMZ FMFHBOU EFTDSJQUJWF
XPSL GPS VT
EP XF DPOTUSVDU UIFTF CBTJT WBSJBCMFT # * EJTQMBZ UIF TJNQMFTU DBTF JO 'ĶĴłĿĲ ƌƉƊ
BQQSPYJNBUF UIF UFNQFSBUVSF EBUB XJUI GPVS EJČFSFOU MJOFBS BQQSPYJNBUJPOT 'JSTU

74. Spline
Basis2
Basis1 Basis3
Basis4
w = [1,–1,1,–1]

75. w = [–1,–1,1,–1]

76. w = [–1,1,1,–1]

77. w = [–1,1,–1,–1]
w = [–1,1,–1,1]

78. Height as a function of age
Obviously not linear
Want function
approximately right
Fit spline as example
But biological model
would do a lot better 0 20 40 60 80
60 80 100 140 180
age (years)
height (cm)

79. Height as a function of age
W
H
S
A
age
weight
height
sex

80. Height as a function of age
W
H
S
A
age
How can age cause anything?
No definable intervention?
Association between age and
height due to accumulated
growth

81. Weighted basis functions
Spline

82. Curves and Splines
Can build very non-linear functions from linear pieces
Splines are powerful geocentric devices
e.g. Average weight only increases with height
e.g. Height increases with age, then levels oﬀ (or declines)
Ideally statistical model has some form as scientific model

83. 0 20 40 60 80
60 80 100 140 180
age (years)
height (cm)

84. 0 20 40 60 80
60 80 100 140 180
age (years)
height (cm)
infancy
childhood
puberty

85. 0 20 40 60 80
60 80 100 140 180
age (years)
height (cm)
infancy
childhood
puberty

86. Course Schedule
Week 1 Bayesian inference Chapters 1, 2, 3
Week 2 Linear models & Causal Inference Chapter 4
Week 3 Causes, Confounds & Colliders Chapters 5 & 6
Week 4 Overfitting / Interactions Chapters 7 & 8
Week 5 MCMC & Generalized Linear Models Chapters 9, 10, 11
Week 6 Integers & Other Monsters Chapters 11 & 12
Week 7 Multilevel models I Chapter 13
Week 8 Multilevel models II Chapter 14
Week 9 Measurement & Missingness Chapter 15
Week 10 Generalized Linear Madness Chapter 16
https://github.com/rmcelreath/stat_rethinking_2023

87. BONUS

88. Full Luxury Bayes
We used two models for two
estimands
But alternative and equivalent
approach is to use one model of
entire causal system
Then use joint posterior to
compute each estimand
W
H
S

89. Full Luxury Bayes
m_SHW_full alist(
# weight
W ~ dnorm(mu,sigma),
mu a[S] ~ dnorm(60,10),
b[S] ~ dunif(0,1),
sigma ~ dunif(0,10),
# height
H ~ dnorm(nu,tau),
nu h[S] ~ dnorm(160,10),
tau ~ dunif(0,10)
), data=dat )
W
H
S

90. Full Luxury Bayes
m_SHW_full alist(
# weight
W ~ dnorm(mu,sigma),
mu a[S] ~ dnorm(60,10),
b[S] ~ dunif(0,1),
sigma ~ dunif(0,10),
# height
H ~ dnorm(nu,tau),
nu h[S] ~ dnorm(160,10),
tau ~ dunif(0,10)
), data=dat )
W
H
S
H
S
W
H
S

91. Full Luxury Bayes
m_SHW_full alist(
# weight
W ~ dnorm(mu,sigma),
mu a[S] ~ dnorm(60,10),
b[S] ~ dunif(0,1),
sigma ~ dunif(0,10),
# height
H ~ dnorm(nu,tau),
nu h[S] ~ dnorm(160,10),
tau ~ dunif(0,10)
), data=dat )
Causal eﬀect is consequence
of intervention
Now simulate each
intervention

92. Simulating Interventions
post Hbar n with( post , {
# simulate W for S=1
H_S1 W_S1 b[,1]*(H_S1-Hbar) , sigma)
# simulate W for S=2
H_S2 W_S2 b[,2]*(H_S2-Hbar) , sigma)
# compute contrast
W_do_S <})
Total causal eﬀect of S on W:
Consequence of changing S at birth
3 4 5 6 7 8 9 10
0.0 0.3 0.6
posterior mean weight contrast (kg)
Density
-20 -10 0 10 20 30
0.00 0.02 0.04
posterior weight contrast (kg)
Density
“p(W|do(S))”

93. Simulating Interventions
post Hbar n with( post , {
# simulate W for S=1
H_S1 W_S1 b[,1]*(H_S1-Hbar) , sigma)
# simulate W for S=2
H_S2 W_S2 b[,2]*(H_S2-Hbar) , sigma)
# compute contrast
W_do_S <})
Total causal eﬀect of S on W:
Consequence of changing S at birth
3 4 5 6 7 8 9 10
0.0 0.3 0.6
posterior mean weight contrast (kg)
Density
-20 -10 0 10 20 30
0.00 0.02 0.04
posterior weight contrast (kg)
Density
“p(W|do(S))”
# automated way
HWsim data=list(S=c(1,2)),
vars=c("H","W"))
W_do_S_auto

94. Inference With Linear Models
With more than two variables, scientific (causal)
model and statistical model not always same
(1) State each estimand
(2) Design unique statistical model for each
(3) Compute each estimand
Or
(1) State each estimand
(2) Compute joint posterior for causal system
(3) Simulate each estimand as an intervention
ONE STAT MODEL
FOR EACH ESTIMAND
ONE SIMULATION
FOR EACH ESTIMAND