Reading Response Discussions

Specify a hypothesis regarding a particular exposure and outcome and a binary effect modifier including specific measures of association (specify the magnitudes of that association you anticipate: I suggest making everything cross-sectional). Using the software of your choice, generate a population with 1000 people under a causal structure consistent with this hypothesis. Draw a simple random sample 100 individuals from this population and estimate the population average exposure-outcome association and the association stratified by your modifier of interest within this subset.  Repeat this 10 times and write the parameter estimates and CI each time.

Exposure (a): rainfall as measured as a deviation from long-term trends, a continuous measure capturing drought/extreme rainfall and normalized

Outcome (y): BMI among children under 5

Effect modifier (w): urban vs rural residence

Data generating rules:

n <- 1000

a <- rnorm(n, 0, 1)

w <- rbinom(n, 1, 0.3)

y_e <- runif(n)

y <-  0.1*a - 0.2*w - 0.1*a*w + y_e

Expected effect sizes:

Among urban (w = 1) = 0

Among rural (w = 0) = 0.1

Prevalence of urban = 0.3

Among total population: 0*0.3 + 0.1*(1-0.3) = 0.07

Results for each sample of 100 individuals

Repeat the data set construction, setting the causal effect to the null.  Again repeat this 10 times and write the parameter estimate and CI each time (if you figure out how to automate it, run it 1000 times and post the histogram of the parameter estimates and p-values).

See attached histograms

Use your code above and also a canned software command to estimate statistical power to detect the difference in means under the settings below:

Note: see attachment for R code and simulation

*n=100, μ₀=.02, μ₁=.12, SD=1, α=.05

Power = 0.07

*n=100,μ₀=.02, μ₁=.12, SD=2, α=.05

Power = 0.04

*n=500, μ₀=.3, μ₁=.3, SD=1, α=.05

Power = 0.025

My apologies but I'm going to skip this week. I attempted to answer but after seeing what others posted I realize I went down the wrong rabbit hole. I will read over everyone's postings to see where I went wrong. Thank you!

Sandeep,

Very nice.  That's a pretty big effect of alcohol and race on SBP in the original simulation (think of it in terms of standard deviations of SBP - does it make sense?)

 Note that although Stata is not thrilled about giving you a power calculation from its canned software you can trick it by moving the two means closer and closer together, i.e., instead of:

. power twomeans 0.3 0.3, sd(1) n(500) alpha(0.05)

use 

. power twomeans 0.3 0.30001, sd(1) n(500) alpha(0.05)

But the point is just to see that even if the two means are identical, there's still an alpha(=0.05) level probability of rejecting the null.

Hi Alice,

Do you mind sharing your code with me for drawing the 1000 samples and storing them in STATA? I did a code for the 1000 samples and I stored it in scalars but I have no idea how to generate an histogram with scalars, maybe I should have used another macro?

Thanks!

Monica

Alice

Nice results - can you post your code?

If you think of power is the probability of rejecting the null under the alternative hypothesis, when the alternative is that the two means are the same, power=alpha.

Usually we think of power as the probability of rejecting the null when it is false assuming the null is that m0=m1 (which is why stata is not happy w/ the power calc when you make m0=m1 the  alternative hypothesis, the same as the null) but sometimes you need to ask questions about power when the null hypothesis is something besides no difference in means.   Suppose for example, that your null hypothesis is that the mean in group 1 is at least 3 points higher than the mean in group 2, and the alternative is that the mean in group 1 is less than 3 points higher than the mean in group 2. You might want to calculate the power to reject the null of a 3+ point difference with any given data set. 

Maria

HW 6

Laura

Yes- I find it a bit easier to simulate continuous outcomes - for binary this is where the log odds idea is useful because you can generate the log odds as a continuous outcome, and then convert the log odds to an odds (by exponentiating) and the odds to a probability (by using prob=odds/(1+odds)). 

As far as choosing an intercept- you can just try some values to see what delivers plausible values.  For continuous outcomes, a mean value (eg 25 for BMI) is a good starting point.  For binary outcomes, it depends how you are generating the variable but if you use the rule I gave you above, a log odds of -3 delivers a probability of about 5%.  A log odds of -2 gives a probability of about 12%. 

Maria

Please see attached.