The questions are in HW6 file, please use R studio only
STA 542: Categorical Data Analysis Ch 7: Loglinear Models (part1) Guangyu Zhu University of Rhode Island Outline 1. Loglinear Model for two-way table 2. Loglinear Models for Three-Way Tables 3. Summary Introduction Log-linear models (or Poisson regression) log(µ) = α+ β1x1 + β2x2 + . . .+ βkxk Interested in patterns of dependence and independence among variables: • Any variables independent? • Strength of associations? • Interactions? 1/50 Loglinear Model for two-way table Example: Belief in Afterlife Rectangle Example: Belief in Afterlife Belief a life after death by race: Race Belief in Afterlife Yes No white 1339 300 black 260 55 other 260 55 2/50 Example: Belief in Afterlife Input contingency table mat <- matrix(c(1339,="" 300,="" 260,="" 55,="" 88,="" 22),="" ncol="2," byrow="T)" dimnames(mat)=""><- list(race="c("white"," "black",="" "other"),="" belief="c("yes"," "no"))="" al.tab=""><- as.table(mat)="" al.tab="" ##="" belief="" ##="" race="" yes="" no="" ##="" white="" 1339="" 300="" ##="" black="" 260="" 55="" ##="" other="" 88="" 22="" function="" odds()="" calculate="" odds="" for="" frequency="" tables="" library(vcd)="" odds(al.tab,="" response="2)" ##="" odds="" for="" race="" by="" belief="" ##="" ##="" white="" black="" other="" ##="" 4.46="" 4.73="" 4.00="" •="" response:="" function="" argument,="" indicating="" the="" response="" variable.="" 3/50="" review:="" pearson's="" chi-squared="" test="" of="" independence="" h0="" :="" x="" and="" y="" independentv.s.="" ha="" :="" x="" and="" y="" dependent="" mles="" under="" h0="" are="" estimated="" expected="" frequencies="" {µ̂ij}.="" µ̂ij="nπ̂i+π̂+j" =="" n="" ni+="" n="" n+j="" n="ni+n+j" n="" usual="" test="" statistic="" is="" pearson's="" chi-squared="" statistic:="" x="" 2="∑" ij="" (nij="" −="" µ̂ij)2="" µ̂ij="∑" all="" cells="" (observed−="" expected)2="" expected="" •="" x="" 2="" has="" a="" large-sample="" chi-squared="" dist.="" under="" h0,="" with="" df="(I" −="" 1)(j="" −="" 1)="" where="" i="number" of="" rows,="" j="number" of="" columns.="" •="" p-value="P(X" 2="" ≥="" x="" 2obs)="right" tail="" probability="" 4/50="" review:="" pearson's="" chi-squared="" test="" of="" independence="" chi-square="" test="" in="" r="" al.chisq=""><- chisq.test(al.tab)="" al.chisq="" ##="" ##="" pearson's="" chi-squared="" test="" ##="" ##="" data:="" al.tab="" ##="" x-squared="0.4," df="2," p-value="0.8" observed="" and="" expected="" frequencies="" al.chisq$observed="" ##="" belief="" ##="" race="" yes="" no="" ##="" white="" 1339="" 300="" ##="" black="" 260="" 55="" ##="" other="" 88="" 22="" al.chisq$expected="" ##="" belief="" ##="" race="" yes="" no="" ##="" white="" 1339.6="" 299.4="" ##="" black="" 257.5="" 57.5="" ##="" other="" 89.9="" 20.1="" 5/50="" loglinear="" model="" for="" two-way="" table="" loglinear="" model="" of="" independence="" motivation="" for="" i="" ×="" j="" table,="" if="" x="" and="" y="" are="" independent="" πij="πi+π+j" for="" all="" i="" ,="" j="" treat="" cell="" counts="" nij="" as="" poisson="" distribution="" with="" mean="" µij="nπij" =="" nπi+π+j="" loglinear="" model="" of="" independence="" log(µij)="λ+" λ="" x="" i="" +="" λ="" y="" j="" •="" λxi="log(πi+):" a="" term="" based="" on="" the="" probability="" in="" row="" i="" •="" λyj="log(π+j" ):="" a="" term="" based="" on="" the="" probability="" in="" column="" j="" .="" 6/50="" loglinear="" model="" of="" independence="" loglinear="" model="" of="" independence="" log(µij)="λ+" λ="" x="" i="" +="" λ="" y="" j="" x="" and="" y="" are="" independent="" ⇔="" loglinear="" model="" of="" independence="" holds="" the="" �tted="" values="" that="" satisfy="" the="" model="" are="" {µ̂ij="ni+n+j/n}" •="" these="" are="" the="" estimated="" expected="" frequencies="" for="" the="" x="" 2="" and="" g="" 2="" tests="" of="" independence.="" •="" those="" tests="" are="" also="" goodness-of-�t="" tests="" of="" this="" loglinear="" model.="" 7/50="" loglinear="" model="" of="" independence:="" interpretation="" suppose="" y="" has="" two="" levels="" (j="2)" in="" the="" loglinear="" model="" of="" independence="" log(µij)="λ+" λ="" x="" i="" +="" λ="" y="" j="" in="" row="" i="" ,="" the="" logit="" for="" p(y="1)" equals="" logit(p(y="1))" =="" log="" [µi1="" µi2="" ]="(λ+" λxi="" +="" λ="" y="" 1="" )−="" (λ+="" λxi="" +="" λy2="" )="λY1" −="" λy2="" 8/50="" loglinear="" model="" of="" independence="" the="" loglinear="" model="" of="" independence="" log(µij)="λ+" λ="" x="" i="" +="" λ="" y="" j="" ⇐⇒="" logit[p(y="1)]" =="" α="" in="" each="" row,="" the="" odds="" of="" response="" in="" column="" 1="" equal="" exp(α)="exp(λY1" −="" λy2="" )="" loglinear="" model="" of="" independence:="" treats="" x="" and="" y="" symmetrically.="" unlike,="" e.g.,="" logistic="" regression="" where="" y="response," x="explanatory." 9/50="" loglinear="" model="" of="" independence="" for="" the="" independence="" model,="" one="" of="" {="" λxi="" }="" is="" redundant,="" and="" one="" of="" {="" λyj="" }="" is="" redundant.="" the="" choice="" of="" constraints="" is="" arbitrary.="" what's="" unique="" about="" the="" parameters?="" •="" the="" di�erences="" between="" them="" are="" unique:(="" λ̂y1="" −="" λ̂y2="" )="unique" value(="" λ̂x1="" −="" λ̂x2="" )="unique" value="" •="" the="" goodness-of-�t="" statistics="" are="" unique="" 10/50="" loglinear="" model="" for="" two-way="" table="" example="" belief="" in="" afterlife:="" loglinear="" model="" example="" belief="" in="" afterlife:="" fitted="" value="" convert="" to="" data.frame="" al=""><- as.data.frame(al.tab)="" al="" ##="" race="" belief="" freq="" ##="" 1="" white="" yes="" 1339="" ##="" 2="" black="" yes="" 260="" ##="" 3="" other="" yes="" 88="" ##="" 4="" white="" no="" 300="" ##="" 5="" black="" no="" 55="" ##="" 6="" other="" no="" 22="" fitting="" independence="" loglinear="" model="" options(contrasts="c("contr.SAS"," "contr.poly"))="" al.mod1=""><- glm(freq="" ~="" race="" +="" belief,="" al,="" family="poisson)" fitted="" value="" of="" independence="" loglinear="" model="" predict(al.mod1,="" type="response" )="" ##="" 1="" 2="" 3="" 4="" 5="" 6="" ##="" 1339.6="" 257.5="" 89.9="" 299.4="" 57.5="" 20.1="" fitted="" value="" in="" chisqured="" test="" al.chisq$expected="" ##="" belief="" ##="" race="" yes="" no="" ##="" white="" 1339.6="" 299.4="" ##="" black="" 257.5="" 57.5="" ##="" other="" 89.9="" 20.1="" 11/50="" example:="" belief="" in="" afterlife�logit="" vs="" loglinear="" loglinear="" model="" coef(al.mod1)="" ##="" (intercept)="" racewhite="" raceblack="" beliefyes="" ##="" 3.00="" 2.70="" 1.05="" 1.50="" logistic="" model="" #="" we="" firstly="" get="" wide="" format="" data="" library(reshape2)="" al_wider=""><- dcast(al,="" race="" ~="" belief,="" value.var="Freq" )="" al_wider="" ##="" race="" yes="" no="" ##="" 1="" white="" 1339="" 300="" ##="" 2="" black="" 260="" 55="" ##="" 3="" other="" 88="" 22="" #="" al.logistic=""><- glm(cbind(yes,="" no)="" ~="" 1,="" family="binomial," al_wider)="" coef(al.logistic)="" ##="" (intercept)="" ##="" 1.5="" so="" the="" odds="" of="" belief="" a="" life="" after="" death="" is="" exp(α)="exp(λ2" −="" λ1)="exp(1.5)" 12/50="" example:="" belief="" in="" afterlife:="" di�erent="" contrasts="" fit="" loglinear="" model="" with="" a="" di�erent="" contrasts="" options(contrasts="c("contr.treatment"," "contr.poly"))="" al.mod2=""><- glm(freq="" ~="" race="" +="" belief,="" al,="" family="poisson)" compare="" two="" sets="" of="" coe�cients="" coef(al.mod1)="" ##="" (intercept)="" racewhite="" raceblack="" beliefyes="" ##="" 3.00="" 2.70="" 1.05="" 1.50="" coef(al.mod2)="" ##="" (intercept)="" raceblack="" raceother="" beliefno="" ##="" 7.20="" -1.65="" -2.70="" -1.50="" black−="" white="1.05−" 2.70="−1.65" fitted="" value="" predict(al.mod1,="" type="response" )="" ##="" 1="" 2="" 3="" 4="" 5="" 6="" ##="" 1339.6="" 257.5="" 89.9="" 299.4="" 57.5="" 20.1="" predict(al.mod2,="" type="response" )="" ##="" 1="" 2="" 3="" 4="" 5="" 6="" ##="" 1339.6="" 257.5="" 89.9="" 299.4="" 57.5="" 20.1="" 13/50="" loglinear="" model="" for="" two-way="" table="" loglinear="" model="" saturated="" model="" loglinear="" model="" saturated="" model="" saturated="" model="" log(µij)="λ+" λ="" x="" i="" +="" λ="" y="" j="" +="" λ="" xy="" ij="" •="" the="" {λxyij="" }="" parameters="" are="" association="" terms="" that="" re�ect="" deviations="" from="" independence.="" •="" the="" independence="" model="" is="" the="" special="" case="" in="" which="" all="" λxyij="0." 14/50="" loglinear="" model:="" number="" of="" parameter="" and="" df="" for="" a="" poisson="" loglinear="" model,="" df="number" of="" poisson="" counts−="" number="" of="" parameters="" for="" an="" i="" ×="" j="" contingency="" table:="" •="" indep.="" model:="" log(µij="" )="λ+" λ="" x="" i="" +="" λ="" y="" j="" no.="" parameters="1+" (i="" −="" 1)="" +="" (j="" −="" 1)="I" +="" j="" −="" 1="" df="IJ" −="" (i="" +="" j="" −="" 1)="(I" −="" 1)(j="" −="" 1)="" •="" saturated="" model:="" log(µij="" )="λ+" λ="" x="" i="" +="" λ="" y="" j="" +="" λ="" xy="" ij="" no.="" parameters="1+" (i="" −="" 1)="" +="" (j="" −="" 1)="" +="" (i="" −="" 1)(j="" −="" 1)="IJ" df="IJ" −="" ij="0" in="" practice,="" unsaturated="" models="" are="" preferred="" 15/50="" example:="" belief="" in="" afterlife="" contingency="" table="" al.tab="" ##="" belief="" ##="" race="" yes="" no="" ##="" white="" 1339="" 300="" ##="" black="" 260="" 55="" ##="" other="" 88="" 22="" fitting="" saturated="" model="" loglinear="" model="" options(contrasts="c("contr.SAS"," "contr.poly"))="" al.satur=""><- glm(freq="" ~="" race="" *="" belief,="" al,="" family="poisson)" coef(al.satur)="" ##="" (intercept)="" racewhite="" raceblack="" beliefyes="" racewhite:beliefyes="" ##="" 3.091="" 2.613="" 0.916="" 1.386="" 0.110="" ##="" raceblack:beliefyes="" ##="" 0.167="" •="" the="" saturated="" model="" gives="" perfect="" �t="" resid(al.satur)="" ##="" [1]="" 0="" 0="" 0="" 0="" 0="" 0="" 16/50="" example:="" belief="" in="" afterlife="" test="" indenpendence="" h0="" :="" independence="" model="" holds="" v.s.="" ha="" :="" saturated="" model="" holds="" lrt="" of="" independence="" anova(al.mod1,="" al.satur,="" test="LRT" )="" ##="" analysis="" of="" deviance="" table="" ##="" ##="" model="" 1:="" freq="" ~="" race="" +="" belief="" ##="" model="" 2:="" freq="" ~="" race="" *="" belief="" ##="" resid.="" df="" resid.="" dev="" df="" deviance="" pr(="">Chi) ## 1 2 0.356 ## 2 0 0.000 2 0.356 0.84 # also report in library(vcd) assocstats(al.tab) ## X^2 df P(> X^2) ## Likelihood Ratio 0.35649 2 0.83674 ## Pearson 0.36008 2 0.83524 ## ## Phi-Coefficient : NA ## Contingency Coeff.: 0.013 ## Cramer's V : 0.013 17/50 Example: Belief in Afterlife�Residuals Pearson residuals ei = yi − µ̂i√ V̂ar(yi ) e <- resid(al.mod1, type = "pearson") e ## 1 2 3 4 5 6 ## -0.0172 0.1581 -0.2012 0.0363 -0.3344 0.4257 for example, for poisson glm the pearson residuals ei = yi − µ̂i√ µ̂i the reason for calling ei a pearson residual is that n∑ i=1 e2i = ∑ i (yi − µ̂i )2 µ̂i = x 2 sum(e^2) ## [1] 0.36 18/50 loglinear model saturated model direct relationships exist between log odds ratios and the {λxyij } association parameters. for example, the model for 2× 2 tables has log odds ratio log θ = log ( µ11µ22 µ12µ21 ) = logµ11 + logµ22 − logµ12 − logµ21 = ( λ+ λx1 + λ y 1 + λ xy 11 ) + ( λ+ λx2 + λ y 2 + λ xy 22 ) − ( λ+ λx1 + λ y 2 + λ xy 12 ) − ( λ+ λx2 + λ y 1 + λ xy 21 ) =λxy11 + λ xy 22 − λxy12 − λxy21 19/50 loglinear model saturated model more generaly: log-odds-ratio comparing levels i and i ′ of x and j and j ′ of y is log(θ) = log (µijµi ′j′ µi ′jµij′ ) = log(µij) + log(µi ′j′)− log(µi ′j)− log(µij′) =(λ+ λxi + λ y j + λ xy ij ) + (λ+ λ x i ′ + λ y j′ + λ xy i ′j′ ) − (λ+ λxi ′ + λyj + λxyi ′j )− (λ+ λxi + λyj′ + λxyij′ ) =λxyij + λ xy i ′j′ − λxyi ′j − λxyij′ substituting the mles of the saturated model (perfect �t) justreproduces the empirical odds ratio nijni ′j′ ni ′jnij′ for the independence model λxyij is 0, and the odds-ratio is e 0 = 1. 20/50 example: belief in afterlife the two interaction terms decide the odds ratio. the estimated odds ratios between belief and race are • exp(0.1096) = 1.12 for white and other, • exp(0.1671) = 1.18 for black and other, • exp(0.1096-0.1671) = 0.94 for white and black. the estimated odds ratios are the same as the sample odds ratios. library(vcd) vcd::oddsratio(al.tab, log = false) ## odds ratios for race and belief ## ## white:black black:other resid(al.mod1,="" type="pearson" )="" e="" ##="" 1="" 2="" 3="" 4="" 5="" 6="" ##="" -0.0172="" 0.1581="" -0.2012="" 0.0363="" -0.3344="" 0.4257="" for="" example,="" for="" poisson="" glm="" the="" pearson="" residuals="" ei="yi" −="" µ̂i√="" µ̂i="" the="" reason="" for="" calling="" ei="" a="" pearson="" residual="" is="" that="" n∑="" i="1" e2i="∑" i="" (yi="" −="" µ̂i="" )2="" µ̂i="X" 2="" sum(e^2)="" ##="" [1]="" 0.36="" 18/50="" loglinear="" model="" saturated="" model="" direct="" relationships="" exist="" between="" log="" odds="" ratios="" and="" the="" {λxyij="" }="" association="" parameters.="" for="" example,="" the="" model="" for="" 2×="" 2="" tables="" has="" log="" odds="" ratio="" log="" θ="log" (="" µ11µ22="" µ12µ21="" )="logµ11" +="" logµ22="" −="" logµ12="" −="" logµ21="(" λ+="" λx1="" +="" λ="" y="" 1="" +="" λ="" xy="" 11="" )="" +="" (="" λ+="" λx2="" +="" λ="" y="" 2="" +="" λ="" xy="" 22="" )="" −="" (="" λ+="" λx1="" +="" λ="" y="" 2="" +="" λ="" xy="" 12="" )="" −="" (="" λ+="" λx2="" +="" λ="" y="" 1="" +="" λ="" xy="" 21="" )="λXY11" +="" λ="" xy="" 22="" −="" λxy12="" −="" λxy21="" 19/50="" loglinear="" model="" saturated="" model="" more="" generaly:="" log-odds-ratio="" comparing="" levels="" i="" and="" i="" ′="" of="" x="" and="" j="" and="" j="" ′="" of="" y="" is="" log(θ)="log" (µijµi="" ′j′="" µi="" ′jµij′="" )="log(µij)" +="" log(µi="" ′j′)−="" log(µi="" ′j)−="" log(µij′)="(λ+" λxi="" +="" λ="" y="" j="" +="" λ="" xy="" ij="" )="" +="" (λ+="" λ="" x="" i="" ′="" +="" λ="" y="" j′="" +="" λ="" xy="" i="" ′j′="" )="" −="" (λ+="" λxi="" ′="" +="" λyj="" +="" λxyi="" ′j="" )−="" (λ+="" λxi="" +="" λyj′="" +="" λxyij′="" )="λXYij" +="" λ="" xy="" i="" ′j′="" −="" λxyi="" ′j="" −="" λxyij′="" substituting="" the="" mles="" of="" the="" saturated="" model="" (perfect="" �t)="" justreproduces="" the="" empirical="" odds="" ratio="" nijni="" ′j′="" ni="" ′jnij′="" for="" the="" independence="" model="" λxyij="" is="" 0,="" and="" the="" odds-ratio="" is="" e="" 0="1." 20/50="" example:="" belief="" in="" afterlife="" the="" two="" interaction="" terms="" decide="" the="" odds="" ratio.="" the="" estimated="" odds="" ratios="" between="" belief="" and="" race="" are="" •="" exp(0.1096)="1.12" for="" white="" and="" other,="" •="" exp(0.1671)="1.18" for="" black="" and="" other,="" •="" exp(0.1096-0.1671)="0.94" for="" white="" and="" black.="" the="" estimated="" odds="" ratios="" are="" the="" same="" as="" the="" sample="" odds="" ratios.="" library(vcd)="" vcd::oddsratio(al.tab,="" log="FALSE)" ##="" odds="" ratios="" for="" race="" and="" belief="" ##="" ##="" white:black="">