[R] Analysing ordinal/nominal data
Piotr Majdak
piotr at majdak.com
Fri Jun 17 12:54:15 CEST 2005
Prof Brian Ripley wrote:
> You have still not defined v' and w', nor the scores (are they estimated
> or not). But the model I suggested is such a model with scores 1,2,...
Sorry for that, here it is:
scores v and w: integer scores, reflecting the ordering of columns/rows.
Agresti suggests to use 1,2,3,... as you did.
v' and w': mean of v and w, respectively
> As I suggested before, a binomial logistic model is appropriate here,
> not a Poisson log-linear one. (They are equivalent, but the binomial
> version is easier to interpret and less wasteful to fit.)
After your suggestions I've read the logit-Chapter in Agresti (1984) and
tried to fit a logit model to my data:
log(m_ij2/m_ij1) = intercept + beta^PR_i(u_i-u') +
beta^ENV_j(v_j-v')
with:
mijk: expected frequencies for cell with indicies i,j,k
u_i: scores of PR with u={0,1,2,3}
v_j: scores of ENV with v={0,1,2}
beta^PR: association parameter for PR with responses
beta^ENV: association parameter for ENV with responses
I wrote in R:
count=c(250,274,285,241,279,299,247,246,255,280,355,362,
230,207,195,239,200,181,233,235,224,200,125,118)
PR=as.integer(gl(4,3,12))-1
ENV=as.integer(gl(3,1,12))-1
countM=matrix(count,12,2)
fit=glm(countM ~ as.integer(PR)+as.integer(ENV),
family=binomial(link=logit))
summary(fit)
R gives me:
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 0.18435 0.07737 2.383 0.017179 *
as.integer(PR) -0.10538 0.03085 -3.416 0.000635 ***
as.integer(ENV) -0.07075 0.04748 -1.490 0.136191
Looking at my data, I know that for PR=0 I have a very weak dependence
on ENV and for PR=3 very strong dependence on ENV. How can I get this
interpretation from the summary above? I'm new in R: is
fit$linearpredictors what I'm looking for?
--
Piotr Majdak
Institut für Schallforschung
Österreichische Akademie der Wissenschaften
Reichsratsstr. 17
A-1010 Wien
Tel.: +43-1-4277-29511
Fax: +43-1-4277-9296
E-Mail: piotr at majdak.com
WWW: http://www.kfs.oeaw.ac.at
More information about the R-help
mailing list