[R] fitting structured conditional (subset) models with loglm

Bert Gunter gunter.berton at gene.com
Mon Feb 10 17:12:09 CET 2014


Perhaps:

?tapply

and/or various wrappers like ?by .

Cheers,
Bert

Bert Gunter
Genentech Nonclinical Biostatistics
(650) 467-7374

"Data is not information. Information is not knowledge. And knowledge
is certainly not wisdom."
H. Gilbert Welch




On Mon, Feb 10, 2014 at 8:05 AM, Michael Friendly <friendly at yorku.ca> wrote:
> With data like the following, a frequency table in data frame form, I'd like
> to fit a collection of loglm models
> of independence of ~ attitude + memory for each combination of education and
> age.
> I can use apply() if I first convert the data to a 2 x 2 x 3 x 3 array, but
> I can't figure out  an
> equivalently simple use of an apply() approach with the data frame form.
>
>> library(MASS)
>> data("Punishment", package = "vcd")
>> str(Punishment)
> 'data.frame':   36 obs. of  5 variables:
>  $ Freq     : num  1 3 20 2 8 4 2 6 1 26 ...
>  $ attitude : Factor w/ 2 levels "no","moderate": 1 1 1 1 1 1 1 1 1 1 ...
>  $ memory   : Factor w/ 2 levels "yes","no": 1 1 1 1 1 1 1 1 1 2 ...
>  $ education: Factor w/ 3 levels "elementary","secondary",..: 1 1 1 2 2 2 3
> 3 3 1 ...
>  $ age      : Factor w/ 3 levels "15-24","25-39",..: 1 2 3 1 2 3 1 2 3 1 ...
>
>> pun <- xtabs(Freq ~ memory + attitude + age + education, data =
>> Punishment)
>>
>> mods.list <- apply(pun, c("age", "education"), function(x) loglm(~memory +
>> attitude, data=x))
>> GSQ <- matrix( sapply(mods.list, function(x)x$lrt), 3, 3)
>> dimnames(GSQ) <- dimnames(mods.list)
>> GSQ
>        education
> age     elementary  secondary       high
>   15-24   4.639061 0.08066111 0.09354563
>   25-39  10.441996 0.96287690 0.48273162
>   40-    12.680802 6.71016542 3.58752829
>> sum(GSQ)
> [1] 39.67937
>
> With the data in data frame format, I can do the same using the subset=
> argument, and
> a series of separate calls (or for loops), but I'd rather us an apply()  (or
> plyr) approach.
>
>> mod.1 <- loglm(Freq ~ memory + attitude, subset=age=="15-24" &
>> education=="elementary", data=Punishment)
>> mod.2 <- loglm(Freq ~ memory + attitude, subset=age=="25-39" &
>> education=="elementary", data=Punishment)
>> mod.3 <- loglm(Freq ~ memory + attitude, subset=age=="40-"   &
>> education=="elementary", data=Punishment)
>> mod.4 <- loglm(Freq ~ memory + attitude, subset=age=="15-24" &
>> education=="secondary", data=Punishment)
>> mod.5 <- loglm(Freq ~ memory + attitude, subset=age=="25-39" &
>> education=="secondary", data=Punishment)
>> mod.6 <- loglm(Freq ~ memory + attitude, subset=age=="40-"   &
>> education=="secondary", data=Punishment)
>> mod.7 <- loglm(Freq ~ memory + attitude, subset=age=="15-24" &
>> education=="high", data=Punishment)
>> mod.8 <- loglm(Freq ~ memory + attitude, subset=age=="25-39" &
>> education=="high", data=Punishment)
>> mod.9 <- loglm(Freq ~ memory + attitude, subset=age=="40-"   &
>> education=="high", data=Punishment)
>>
>> mod.list <- list(mod.1, mod.2,mod.3, mod.4, mod.5, mod.6, mod.7, mod.8,
>> mod.9)
>>
>> GSQ <- matrix( sapply(mod.list, function(x)x$lrt), 3, 3)
>> dimnames(GSQ) <- list(age = levels(Punishment$age),
> +                       education = levels(Punishment$education)
> +                       )
>> GSQ
>        education
> age     elementary  secondary       high
>   15-24   4.639061 0.08066111 0.09354563
>   25-39  10.441996 0.96287690 0.48273162
>   40-    12.680802 6.71016542 3.58752829
>>
>
> --
> Michael Friendly     Email: friendly AT yorku DOT ca
> Professor, Psychology Dept. & Chair, Quantitative Methods
> York University      Voice: 416 736-2100 x66249 Fax: 416 736-5814
> 4700 Keele Street    Web:http://www.datavis.ca
> Toronto, ONT  M3J 1P3 CANADA
>
> ______________________________________________
> 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