Adaikalavan Ramasamy
2010-Jul-13 18:32 UTC
[R] working out main effect variance when different parameterization is used and interaction term exists
Dear all, Apologies if this question is bit theoretical and for the longish email. I am meta-analyzing the coefficients and standard errors from multiple studies where the raw data is not available. Each study analyst runs a model that includes an interaction term for, say, between sex and smoking and age. Here is an illustrative example example for one study: set.seed(1066) status <- rbinom( 1000, 1, 0.2 ) males <- rbinom( 1000, 1, 0.6 ) smoke <- rbinom( 1000, 1, 0.3 ) age <- runif(1000, min=20, max=80) coef( summary( f1 <- glm( status ~ males*smoke + age, family="binomial" ) ) ) # Estimate Std. Error z value Pr(>|z|) # (Intercept) -1.520399871 0.284464584 -5.3447774 9.052825e-08 # males 0.213851446 0.201717381 1.0601538 2.890746e-01 # smoke -0.123103049 0.292346483 -0.4210861 6.736922e-01 # age -0.001056007 0.004612947 -0.2289223 8.189293e-01 # males:smoke 0.283775173 0.362821438 0.7821345 4.341355e-01 Now, unfortunately some analysts coded sex as females instead of males. Using the same dataset, I get the following output with females: females <- 1 - males coef( summary( f1 <- glm( status ~ females*smoke + age, family="binomial" )) ) # Estimate Std. Error z value Pr(>|z|) # (Intercept) -1.306548425 0.262573162* -4.9759405 6.493160e-07 # females -0.213851446 0.201717381* -1.0601538 2.890746e-01 # smoke 0.160672124 0.214923130* 0.7475795 4.547138e-01 # age -0.001056007 0.004612947 -0.2289223 8.189293e-01 # females:smoke -0.283775173 0.362821438 -0.7821345 4.341355e-01 I have worked out algebrically (and numerically) the following: Beta(females) = -Beta(males) Var(females) = Var(males) Beta(females:smoke) = -Beta(males:smoke) Var(females:smoke) = Var(males:smoke) Beta(smoke | fit1) = Beta(smoke | fit2) + Beta(females:smoke) = 0.160672124 -0.283775173 = -0.1231030 How can I calculate the Var(smoke | fit1) from Var(smoke | fit2) ? I tried to derive this algebrically but ended up with a covariance term which I could not solve. If I could cleverly convert Var(smoke | fit2) to Var(smoke | fit1) then I could avoid going back to each analyst since this particular analyses is only one of many hundreds we run and it would be annoying for each analyst to use the same parameterisation. Any suggestions is much appreciated. Many thanks in advance. Regards, Adai