[R] new to R need urgent help!

Abhijit Dasgupta adasgupta at araastat.com
Fri Jun 24 03:02:59 CEST 2011


On Jun 23, 2011, at 4:42 PM, elisheva corn <elishevacorn at gmail.com> wrote:

> hi all-
> 
> I am doing some research, have never used R before until today and need to
> understand the following program for a project.
> if some one could PLEASE help me understand this program ASAP i would
> GREATLY appreciate it (any syntax/ statistic comments would be great)
> 
> PLEASE PLEASE HELP!!  THANKYOU!!!
> -on a side note, it seems to me that R doesnt include the pv, and it was
> calculated seperatly, is this true?
> 
> 
> fit=gee(foci~as.factor(time)*cond,id=exper,data=drt,family=poisson(link =
> "log"))
You apparently have count data (foci) which is measured repeatedly within exper, and you're interested in how foci changes with time and condition including their interaction. The code fits a generalized estimating equation (GEE) model, which can be an appropriate model for repeated measures data. See, for example, Diggle, Liang, Zeger & Heagerty for background. 
> Beginning Cgee S-function, @(#) geeformula.q 4.13 98/01/27
> running glm to get initial regression estimate
>                  (Intercept)             as.factor(time)24
>                     3.051177                     -2.705675
>                  condHypoxia as.factor(time)24:condHypoxia
>                    -0.402259                      1.429034
>> pv=2*(1-pnorm(abs(summary(fit)$coef[,5])))
>> data.frame(summary(fit)$coef,pv)

The gee package doesn't compute the value directly, though other functions like lm, glm and others do.  What the code does is use the robust z statistic, which is the estimate/robust se, and relate it to the standard normal distribution. 
>                               Estimate Naive.S.E.   Naive.z Robust.S.E.
> Robust.z
> (Intercept)                    3.051177 0.02221052 137.37527  0.04897055
> 62.306363
> as.factor(time)24             -2.705675 0.10890056 -24.84537  0.19987174
> -13.537057
> condHypoxia                   -0.402259 0.03907961 -10.29332  0.10661248
> -3.773095
> as.factor(time)24:condHypoxia  1.429034 0.12549576  11.38711  0.17867421
> 7.997988
>                                        pv
> (Intercept)                   0.000000e+00
> as.factor(time)24             0.000000e+00
> condHypoxia                   1.612350e-04
> as.factor(time)24:condHypoxia 1.332268e-15
>> ftable(table(drt$cond,drt$time,predict(fit)))
>             0.345501643340608 1.37227675004058 2.64891772174934
> 3.05117673373261
> 
> 
> Oxia    0.5                  0                0
> 0              485
>        24                 315                0
> 0                0
> Hypoxia 0.5                  0                0
> 346                0
>        24                   0              449
> 0                0
>> ## 3-th term gives the difference between the Hypoxia/Oxia at time=0.5
>> ## the difference between Hypoxia/Oxia at time=24
>> L=matrix(c(0,0,1,1),nrow=1)
>> fit$coef[L==1]
>                  condHypoxia as.factor(time)24:condHypoxia
>                    -0.402259                      1.429034
>> L%*%fit$coef
>         [,1]
> [1,] 1.026775
>> wald.test(fit$robust.variance,fit$coef,L=L)
> Wald test:
> ----------
> 
> Chi-squared test:
> X2 = 23.8, df = 1, P(> X2) = 1.1e-06
>> 
> 
>    [[alternative HTML version deleted]]
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.



More information about the R-help mailing list