# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: Python 3 (ipykernel)
#     language: python
#     name: python3
# ---

# ## Problem set #2
#
#
# ### Q1
# ##### Which of the linear models below are multiple linear regressions?
#
# (1) µi = α + βxi
#
# (2) µi = β1 * xi + β2 * zi
#
# (3) µi = α + β(xi − zi)
#
# (4) µi = α + β1 * xi + β2 * zi

# cd C:\Users\edsea\My Drive\Backup\UM6P\Data Science 2 class\DatSci2_2022_repo

import matplotlib.pyplot as plt
import scipy.stats as st
import numpy as np
import pandas as pd
import pymc3 as pm
from quap import quap
import arviz as az
import statsmodels.api as sm
import math
d = pd.read_csv("./Data/drinking_fatalities.csv", sep=";", header=0)

# ### Q2
# I am trying to model the effect of education on fame - in other words I want to know whether having more years of education improves a person's chances of being famous.
#
# Give 1 example of a variable I might want to control for in my model. What could happen if I don't control for it?

# ### Q3
# Write the model definition for the above model. Should you center any variables? Justify your priors.

# ### Q4
# I want to set up an experiment looking at the effect of social media use on anxiety and test scores. Help me simulate this.
#
# First simulate baseline anxiety levels for a group of 100 students. Anxiety is scored between 0 and 10 and the average score is 5 with a standard deviation of 2. Remember that you can't score lower than 0 or higher than 10!
#
# Next assign half of them at random to the experimental group which is given access to social media. 
#
# Simulate student's scores after the treatment. On average, the control group (with no access to social media) have the same score, while the experimental group increased their anxiety by 1 point on average.
#
# Simulate student's test scores out of 100. The average score for a student without anxiety is 80% with a standard deviation of of 15 points, but each point of anxiety decreases a student's score by 2 points on average. (again, remember that scores cannot be more than 100 or less than 0).
#
# Create a data frame with experimental status, post-experiment anxiety, and test scores and show the first 6 rows.

# +
N =100
a0 = np.random.normal(5,2,size=N)
a0[a0<0] = 0
a0[a0>10] = 10
SM = np.repeat([0,1],N/2)
a1 = np.random.normal(a0 + SM*1,0.1)
a1[a1<0]=0
a1[a1>10]=10
TS = np.random.normal(80-a1*2, 15)
TS[TS<0] = 0
TS[TS>100] = 100

df = pd.DataFrame.from_dict({"Social_Media" : SM,"Anxiety" : a1, "Test_score" : TS})
df.head()
# -

# ### Q5
# Your colleague tells you to put model test scores using a multiple linear regression with anxiety and experimental status as predictor variables. Why is this a bad idea?
#
# To appease them, you do it anyway. Present the results in a forest plot. Interpret the meaning of each of the parameter estimates. What model should you run instead? Present these results in another forest plot and interpret the results for your colleague.

# +
with pm.Model() as anxiety_model:
    
    B1 = pm.Normal("B1",0,1)
    B2 = pm.Normal("B2",0,1)
    a = pm.Normal("a",50,20)
    sigma = pm.Exponential("sigma",1)
    mu = a + B1*Anxiety + B2*Social_Media
    
    Test_score = np.Normal("Test_score",mu,sigma, observed=df.Test_score)
    data,dist = quap(vars=[a,B1,B2,sigma],n_samples=10_000)
    
    
# -

# ### Q6
# Load the data /Data/drinking_fatalities.csv. This is a dataset of US states between 1982 and 1988, with data on the minimum legal drinking age. Your colleague from the previous question tells you they have found an interesting result. They find that states with more average miles driven per person have higher minimum legal drinking ages. They interpret this to mean that states adjust their legal drinking ages upwards when people drive a lot.
#
# Run a simple linear model with "mlda" as the outcome and "avgmiles" as the predictor. Present the parameter estimates and their 89% hdis and interpret each of them. Plot the predicted relationship between avgmiles and mlda along with the 89% hdi for the mean and the 89% hdi for the predictions.
#
# Do you come to the same conclusion as your colleague?

# ### Q7
# You ask your colleague how they specified their model. They tell you they controlled for the traffic fatality rate, trafficmort. Rerun the previous model adding trafficmort as a control. Once again, present the parameter estimates and interpret them, and plot the predicted relationship between avgmiles and mlda along with the 89% hdi for the mean and the 89% hdi for the predictions. What is different in this model?
#
# Think about the causal relationship between traffic fatality rates, average driving rates, and legal drinking age. Do you think your colleague's model specification was appropriate? Why or why not? What do you conclude about the true relationship between driving rates and the drinking age?



StarcoderdataJupyterScriptsDedupFiltered