Because the ability estimates are based on different functions of the data. With Bigsteps, the ability estimate (theta) is found by maximizing the data likelihood wrt theta, call it f(x). But, in the MAP method of estimation, ability estimates are derived from f(x)g(x) where g(x) ~ N(\mu, \sigma^2).

There is technical documentation distributed with MiscPsycho that shows what the likelihoods are in MAP and ML estimation. In the former, you maximize the posterior and in the latter you maximize only the data likelihood. 


-----Original Message-----
From: Iasonas Lamprianou [mailto:lamprianou@yahoo.com]
Sent: Sun 2/22/2009 11:49 AM
To: Douglas Bates; Doran, Harold
Cc: r-sig-mixed-models@r-project.org
Subject: RE: [R-sig-ME] Rasch model
 
Thank you DH,
I used bigsteps and Analysis which both give exactly the same results. The ltm package also gave similar results to the BigSteps and Analysis. However, why does lme4 give different ability estimates, while it gives practically the same difficulty estimates?
In case you are interested, I attach the raw data as well as a comparison with the results of various packages in an excel file

thanks

Dr. Iasonas Lamprianou
Department of Education
The University of Manchester
Oxford Road, Manchester M13 9PL, UK
Tel. 0044  161 275 3485
iasonas.lamprianou@manchester.ac.uk


--- On Sun, 22/2/09, Doran, Harold <HDoran@air.org> wrote:

> From: Doran, Harold <HDoran@air.org>
> Subject: RE: [R-sig-ME] Rasch model
> To: "Douglas Bates" <bates@stat.wisc.edu>, lamprianou@yahoo.com
> Cc: r-sig-mixed-models@r-project.org
> Date: Sunday, 22 February, 2009, 2:33 PM
> What "traditional" method (and what packages) did
> you use for ability estimation in the other packages? When
> you treat students as random as you did and then get the
> ability estimates, this is equivalent to MAP estimation.
> Most traditional packages use ML estimation, which limits
> one's ability to get estimates for those individuals
> with perfect or 0 scores since the likelihood function is
> unbounded. Arbitrary scoring rules are typically applied in
> these cases.
> 
> You might check out the irt.ability() function in the
> MiscPsycho package for generating abilities as it offer you
> MAP, EAP and ML estimates in R.
> 
> 
> -----Original Message-----
> From: r-sig-mixed-models-bounces@r-project.org on behalf of
> Douglas Bates
> Sent: Sat 2/21/2009 9:19 AM
> To: lamprianou@yahoo.com
> Cc: r-sig-mixed-models@r-project.org
> Subject: Re: [R-sig-ME] Rasch model
>  
> On Sat, Feb 21, 2009 at 5:31 AM, Iasonas Lamprianou
> <lamprianou@yahoo.com> wrote:
> > Dear friends,
> > I am running a simple model where 424 pupils took 13
> test items. I model the item easiness as fixed effects and
> the pupil abilities as random effects (this is basically the
> 'so-called' marginal maximum likelihood estimation
> used in some Rasch software packages. I also run the
> analysis using two 'traditional' Rasch packages. I
> have found that there is a near-perfect correlation between
> item estimates from the 'traditional' packages and
> lme4. See the results below. However, when I use
> ranef(model$id) to  get the 'ability' estimates of
> the pupils, the correlation is just 0.9! Shouldnt the
> correlation be much bigger? I mean, how would you estimate
> the ability estimates of the pupils in this context?
> 
> > Thanks for any help
> 
> As far as I know the conditional modes of the random
> effects for
> students should correspond to the student ability estimates
> from
> marginal maximum likelihood.  However, I don't know
> much about
> commercial Rasch packages calculate these parameters.
> 
> With more items is becomes reasonable to model both the
> item easiness
> and the pupil abilities as random effects which is not
> easily done in
> commercial packages.  There is a chapter in the book
> "Exploratory Item
> Response Models" describing the model and why it is
> desirable but
> without indication of how it could be fit.  In the paper
> 
> @article{Dowling:Bliese:Bates:Doran:2007:JSSOBK:v20i02,
>   author =	"Harold Doran and Douglas Bates and Paul
> Bliese and Maritza
>  Dowling",
>   title =	"Estimating the Multilevel Rasch Model: With
> the lme4 Package",
>   journal =	"Journal of Statistical Software",
>   volume =	"20",
>   number =	"2",
>   pages =	"1--18",
>   day =  	"22",
>   month =	"2",
>   year = 	"2007",
>   CODEN =	"JSSOBK",
>   ISSN = 	"1548-7660",
>   bibdate =	"2007-02-22",
>   URL =  	"http://www.jstatsoft.org/v20/i02",
>   accepted =	"2007-02-22",
>   acknowledgement = "",
>   keywords =	"",
>   submitted =	"2006-10-01",
> }
> 
> we show how to fit that model and generalizations of it
> where items
> and/or subjects fall into groups.  Also, there are several
> slides in
> the section "Item Response Models as GLMMs" of
> the workshop
> presentation at
> http://www.stat.wisc.edu/~bates/UseR2008/WorkshopD.pdf
> related to fitting Rasch models with crossed random
> effects.
> 
> 
> 
> > Generalized linear mixed model fit by the Laplace
> approximation
> > Formula: score ~ 0 + item + (1 | id)
> >   Data: rasch_data
> >  AIC  BIC logLik deviance
> >  6502 6595  -3237     6474
> > Random effects:
> >  Groups Name        Variance Std.Dev.
> >  id     (Intercept) 1.9821   1.4079
> > Number of obs: 5512, groups: id, 424
> >
> > Fixed effects:
> >         Estimate Std. Error z value Pr(>|z|)
> > item   1  0.05819    0.13004   0.447 0.654524
> > item   2  1.26791    0.14126   8.976  < 2e-16 ***
> > item   3  0.93972    0.13615   6.902 5.12e-12 ***
> > item   4  0.70219    0.13345   5.262 1.43e-07 ***
> > item   5  0.36877    0.13099   2.815 0.004874 **
> > item   6  0.52699    0.13197   3.993 6.52e-05 ***
> > item   7 -0.79568    0.13404  -5.936 2.92e-09 ***
> > item   8  0.38186    0.13106   2.914 0.003572 **
> > item   9 -0.61867    0.13239  -4.673 2.97e-06 ***
> > item  10  0.30358    0.13069   2.323 0.020184 *
> > item  11  0.23870    0.13044   1.830 0.067262 .
> > item  12 -0.79568    0.13404  -5.936 2.92e-09 ***
> > item  13 -0.47278    0.13137  -3.599 0.000320 ***
> > ---
> > Signif. codes:  0 '***' 0.001 '**'
> 0.01 '*' 0.05 '.' 0.1 ' ' 1
> >
> > Correlation of Fixed Effects:
> >         item1 item2 item3 item4 item5 item6 item7
> item8 item9 item10 item11 item12
> > item   2 0.258
> > item   3 0.269 0.252
> > item   4 0.275 0.256 0.265
> > item   5 0.281 0.258 0.268 0.274
> > item   6 0.278 0.257 0.267 0.273 0.277
> > item   7 0.273 0.243 0.255 0.262 0.269 0.266
> > item   8 0.280 0.258 0.268 0.274 0.279 0.277 0.269
> > item   9 0.277 0.248 0.259 0.266 0.273 0.270 0.271
> 0.273
> > item  10 0.281 0.258 0.269 0.275 0.280 0.278 0.270
> 0.280 0.274
> > item  11 0.282 0.258 0.269 0.275 0.280 0.278 0.271
> 0.280 0.275 0.281
> > item  12 0.273 0.243 0.255 0.262 0.269 0.266 0.268
> 0.269 0.271 0.270  0.271
> > item  13 0.279 0.251 0.263 0.269 0.276 0.273 0.272
> 0.276 0.275 0.277  0.278  0.272
> >> MML_Rasch <- lmer(score ~ 0+item+(1|id),
> rasch_data, family=binomial(link="logit"))
> >> summary(MML_Rasch)
> > Generalized linear mixed model fit by the Laplace
> approximation
> > Formula: score ~ 0 + item + (1 | id)
> >   Data: rasch_data
> >  AIC  BIC logLik deviance
> >  6502 6595  -3237     6474
> > Random effects:
> >  Groups Name        Variance Std.Dev.
> >  id     (Intercept) 1.9821   1.4079
> > Number of obs: 5512, groups: id, 424
> >
> > Fixed effects:
> >         Estimate Std. Error z value Pr(>|z|)
> > item   1  0.05819    0.13004   0.447 0.654524
> > item   2  1.26791    0.14126   8.976  < 2e-16 ***
> > item   3  0.93972    0.13615   6.902 5.12e-12 ***
> > item   4  0.70219    0.13345   5.262 1.43e-07 ***
> > item   5  0.36877    0.13099   2.815 0.004874 **
> > item   6  0.52699    0.13197   3.993 6.52e-05 ***
> > item   7 -0.79568    0.13404  -5.936 2.92e-09 ***
> > item   8  0.38186    0.13106   2.914 0.003572 **
> > item   9 -0.61867    0.13239  -4.673 2.97e-06 ***
> > item  10  0.30358    0.13069   2.323 0.020184 *
> > item  11  0.23870    0.13044   1.830 0.067262 .
> > item  12 -0.79568    0.13404  -5.936 2.92e-09 ***
> > item  13 -0.47278    0.13137  -3.599 0.000320 ***
> > ---
> > Signif. codes:  0 '***' 0.001 '**'
> 0.01 '*' 0.05 '.' 0.1 ' ' 1
> >
> > Correlation of Fixed Effects:
> >         item1 item2 item3 item4 item5 item6 item7
> item8 item9 item10 item11 item12
> > item   2 0.258
> > item   3 0.269 0.252
> > item   4 0.275 0.256 0.265
> > item   5 0.281 0.258 0.268 0.274
> > item   6 0.278 0.257 0.267 0.273 0.277
> > item   7 0.273 0.243 0.255 0.262 0.269 0.266
> > item   8 0.280 0.258 0.268 0.274 0.279 0.277 0.269
> > item   9 0.277 0.248 0.259 0.266 0.273 0.270 0.271
> 0.273
> > item  10 0.281 0.258 0.269 0.275 0.280 0.278 0.270
> 0.280 0.274
> > item  11 0.282 0.258 0.269 0.275 0.280 0.278 0.271
> 0.280 0.275 0.281
> > item  12 0.273 0.243 0.255 0.262 0.269 0.266 0.268
> 0.269 0.271 0.270  0.271
> > item  13 0.279 0.251 0.263 0.269 0.276 0.273 0.272
> 0.276 0.275 0.277  0.278  0.272
> >
> >
> > Dr. Iasonas Lamprianou
> > Department of Education
> > The University of Manchester
> > Oxford Road, Manchester M13 9PL, UK
> > Tel. 0044  161 275 3485
> > iasonas.lamprianou@manchester.ac.uk
> >
> >
> > --- On Sat, 21/2/09,
> r-sig-mixed-models-request@r-project.org
> <r-sig-mixed-models-request@r-project.org> wrote:
> >
> >> From: r-sig-mixed-models-request@r-project.org
> <r-sig-mixed-models-request@r-project.org>
> >> Subject: R-sig-mixed-models Digest, Vol 26, Issue
> 33
> >> To: r-sig-mixed-models@r-project.org
> >> Date: Saturday, 21 February, 2009, 11:00 AM
> >> Send R-sig-mixed-models mailing list submissions
> to
> >>       r-sig-mixed-models@r-project.org
> >>
> >> To subscribe or unsubscribe via the World Wide
> Web, visit
> >>      
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> >> or, via email, send a message with subject or body
> >> 'help' to
> >>       r-sig-mixed-models-request@r-project.org
> >>
> >> You can reach the person managing the list at
> >>       r-sig-mixed-models-owner@r-project.org
> >>
> >> When replying, please edit your Subject line so it
> is more
> >> specific
> >> than "Re: Contents of R-sig-mixed-models
> >> digest..."
> >>
> >>
> >> Today's Topics:
> >>
> >>    1. Re: Power calculations for random effect
> models (Greg
> >> Snow)
> >>
> >>
> >>
> ----------------------------------------------------------------------
> >>
> >> Message: 1
> >> Date: Fri, 20 Feb 2009 18:12:52 -0700
> >> From: Greg Snow <Greg.Snow@imail.org>
> >> Subject: Re: [R-sig-ME] Power calculations for
> random
> >> effect models
> >> To: "christos@nuverabio.com"
> >> <christos@nuverabio.com>,
> >>       "r-sig-mixed-models@r-project.org"
> >> <r-sig-mixed-models@r-project.org>
> >> Message-ID:
> >>      
> <B37C0A15B8FB3C468B5BC7EBC7DA14CC61CA3E2EBE@LP-EXMBVS10.CO.IHC.COM>
> >> Content-Type: text/plain;
> charset="us-ascii"
> >>
> >> Since you already have a function to simulate the
> data, it
> >> should be fairly simple to use simulations to
> compare the
> >> methods you suggest below, or to do a simple test
> based on
> >> the likelihood ratio, but simulating the cutoff
> rather than
> >> depending on the chi-squared approximation.
> >>
> >> See:
> >>
> http://finzi.psych.upenn.edu/R/Rhelp08/archive/156522.html
> >>
> >> for some examples of assessing these types of
> questions
> >> using simulation (I have an updated version of the
> last
> >> example if you want it as well).
> >>
> >> Hope this helps,
> >>
> >> --
> >> Gregory (Greg) L. Snow Ph.D.
> >> Statistical Data Center
> >> Intermountain Healthcare
> >> greg.snow@imail.org
> >> 801.408.8111
> >>
> >>
> >> > -----Original Message-----
> >> > From:
> r-sig-mixed-models-bounces@r-project.org
> >> [mailto:r-sig-mixed-
> >> > models-bounces@r-project.org] On Behalf Of
> Christos
> >> Hatzis
> >> > Sent: Friday, February 20, 2009 11:01 AM
> >> > To: r-sig-mixed-models@r-project.org
> >> > Subject: [R-sig-ME] Power calculations for
> random
> >> effect models
> >> >
> >> > Hello,
> >> >
> >> > I have a nested random effects model of the
> form
> >> >
> >> > Yijk = M + Ai + Bj(i) + Ek(ij)
> >> >
> >> > where biopsies (B) are nested within persons
> (A) and
> >> arrays (E) are
> >> > nested
> >> > within biopsies.
> >> > I am interested in estimating the power of a
> study of
> >> a given size to
> >> > determine whether the variance associated
> with B is
> >> trivial, i.e. H0:
> >> > var_b
> >> > = 0 vs Ha: var_b > 0 at a fixed Type-I
> error rate.
> >> >
> >> > I have written a function to simulate data
> from this
> >> model and used
> >> > lmer to
> >> > estimate the random effects as shown below. 
> What is
> >> the recommended
> >> > way of
> >> > going about testing the null hypothesis?  One
> approach
> >> would be to use
> >> > the
> >> > estimated standard deviation of the random
> effect and
> >> assuming
> >> > normality to
> >> > test whether the appropriate CI contains
> zero. Another
> >> would be to use
> >> > the
> >> > theoretical chi-square distribution for
> var_b, but
> >> that would require
> >> > the
> >> > appropriate degrees of freedom.  Or to use
> mcmc to
> >> estimate the
> >> > distribution
> >> > of var_b and use this distribution for
> inference.
> >> >
> >> > I would think that the mcmc approach is the
> >> recommended one, but I
> >> > would
> >> > appreciate any advise on this. In this case,
> I have
> >> tried to run the
> >> > MCMC
> >> > simulation (which runs fine), but have not
> been able
> >> yet to figure out
> >> > how
> >> > to use the results of MCMC to test the above
> >> hypothesis.  Any hints or
> >> > pointing out materials that explain how to
> use MCMC
> >> for inference on
> >> > random
> >> > effects will be very much appreciated.
> >> >
> >> > Thank you.
> >> > -Christos
> >> >
> >> > Christos Hatzis, Ph.D.
> >> > Nuvera Biosciences, Inc.
> >> > 400 West Cummings Park
> >> > Suite 5350
> >> > Woburn, MA 01801
> >> > Tel: 781-938-3830
> >> > www.nuverabio.com
> <http://www.nuverabio.com/>
> >> >
> >> > > fake.dt <- tumor.heter.dt(10, 3, 2,
> k=0, sa=4,
> >> sb=.6, se=.2)
> >> > > fake.dt
> >> >    person biopsy array      resp bioWper
> >> > 1       1      1     1  4.308434     1:1
> >> > 2       1      1     2  4.293186     1:1
> >> > 3       1      2     1  5.503841     1:2
> >> > 4       1      2     2  5.640362     1:2
> >> > 5       1      3     1  5.579461     1:3
> >> > 6       1      3     2  5.416201     1:3
> >> > 7       2      1     1 12.479513     2:1
> >> > 8       2      1     2 12.311426     2:1
> >> > 9       2      2     1 13.283566     2:2
> >> > 10      2      2     2 13.138276     2:2
> >> > 11      2      3     1 12.954277     2:3
> >> > 12      2      3     2 13.059925     2:3
> >> > 13      3      1     1 11.649726     3:1
> >> > 14      3      1     2 11.694472     3:1
> >> > 15      3      2     1 12.342050     3:2
> >> > 16      3      2     2 12.316214     3:2
> >> > 17      3      3     1 12.762695     3:3
> >> > 18      3      3     2 12.451338     3:3
> >> > 19      4      1     1 17.168248     4:1
> >> > 20      4      1     2 17.573315     4:1
> >> > 21      4      2     1 16.885176     4:2
> >> > 22      4      2     2 16.819536     4:2
> >> > 23      4      3     1 15.924120     4:3
> >> > 24      4      3     2 16.038491     4:3
> >> > 25      5      1     1  7.981279     5:1
> >> > 26      5      1     2  7.777929     5:1
> >> > 27      5      2     1  6.701801     5:2
> >> > 28      5      2     2  6.978284     5:2
> >> > 29      5      3     1  7.136417     5:3
> >> > 30      5      3     2  7.378498     5:3
> >> > 31      6      1     1 21.499796     6:1
> >> > 32      6      1     2 21.397173     6:1
> >> > 33      6      2     1 22.551116     6:2
> >> > 34      6      2     2 22.521006     6:2
> >> > 35      6      3     1 21.027998     6:3
> >> > 36      6      3     2 21.416872     6:3
> >> > 37      7      1     1 13.312857     7:1
> >> > 38      7      1     2 13.559062     7:1
> >> > 39      7      2     1 13.753016     7:2
> >> > 40      7      2     2 13.608642     7:2
> >> > 41      7      3     1 13.556446     7:3
> >> > 42      7      3     2 13.400672     7:3
> >> > 43      8      1     1 12.758548     8:1
> >> > 44      8      1     2 12.486574     8:1
> >> > 45      8      2     1 13.388409     8:2
> >> > 46      8      2     2 13.263029     8:2
> >> > 47      8      3     1 12.991308     8:3
> >> > 48      8      3     2 12.962116     8:3
> >> > 49      9      1     1 11.420214     9:1
> >> > 50      9      1     2 11.308010     9:1
> >> > 51      9      2     1 13.186774     9:2
> >> > 52      9      2     2 12.778966     9:2
> >> > 53      9      3     1 11.625079     9:3
> >> > 54      9      3     2 11.637015     9:3
> >> > 55     10      1     1  9.108635    10:1
> >> > 56     10      1     2  8.895658    10:1
> >> > 57     10      2     1  9.718046    10:2
> >> > 58     10      2     2  9.552510    10:2
> >> > 59     10      3     1  8.794893    10:3
> >> > 60     10      3     2  9.047379    10:3
> >> > > heter.lmer <- lmer(resp ~ (1 |
> person) + (1 |
> >> bioWper), fake.dt)
> >> > > heter.lmer
> >> > Linear mixed model fit by REML
> >> > Formula: resp ~ (1 | person) + (1 | bioWper)
> >> >    Data: fake.dt
> >> >    AIC   BIC logLik deviance REMLdev
> >> >  98.24 106.6 -45.12    92.85   90.24
> >> > Random effects:
> >> >  Groups   Name        Variance  Std.Dev.
> >> >  bioWper  (Intercept)  0.312327 0.55886
> >> >  person   (Intercept) 21.780993 4.66701
> >> >  Residual              0.020633 0.14364
> >> > Number of obs: 60, groups: bioWper, 30;
> person, 10
> >> >
> >> > Fixed effects:
> >> >             Estimate Std. Error t value
> >> > (Intercept)   12.368      1.479    8.36
> >> >
> >> >
> >> > >
> VarCorr(heter.lmer)[["bioWper"]]
> >> >             (Intercept)
> >> > (Intercept)   0.3123273
> >> > attr(,"stddev")
> >> > (Intercept)
> >> >   0.5588625
> >> > attr(,"correlation")
> >> >             (Intercept)
> >> > (Intercept)           1
> >> >
> >> > > heter.mcmc <- mcmcsamp(heter.lmer,
> n=1000)
> >> > > str(heter.mcmc)
> >> > Formal class 'merMCMC' [package
> >> "lme4"] with 9 slots
> >> >   ..@ Gp      : int [1:3] 0 30 40
> >> >   ..@ ST      : num [1:2, 1:1000] 3.89 32.49
> 1.44
> >> 27.06 1.24 ...
> >> >   ..@ call    : language lmer(formula = resp
> ~ (1 |
> >> person) + (1 |
> >> > bioWper),
> >> > data = fake.dt)
> >> >   ..@ deviance: num [1:1000] 92.9 92.9 114.5
> 119.1 134
> >> ...
> >> >   ..@ dims    : Named int [1:18] 2 60 1 40 1
> 2 0 1 2 5
> >> ...
> >> >   .. ..- attr(*, "names")= chr
> [1:18]
> >> "nt" "n" "p"
> "q" ...
> >> >   ..@ fixef   : num [1, 1:1000] 12.4 14.1
> 11.9 12.6
> >> 12.6 ...
> >> >   .. ..- attr(*, "dimnames")=List
> of 2
> >> >   .. .. ..$ : chr "(Intercept)"
> >> >   .. .. ..$ : NULL
> >> >   ..@ nc      : int [1:2] 1 1
> >> >   ..@ ranef   : num[1:40, 0 ]
> >> >   ..@ sigma   : num [1, 1:1000] 0.144 0.128
> 0.126
> >> 0.212 0.228 ...
> >> >
> >> >
> _______________________________________________
> >> > R-sig-mixed-models@r-project.org mailing list
> >> >
> >>
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> >>
> >>
> >>
> >> ------------------------------
> >>
> >> _______________________________________________
> >> R-sig-mixed-models mailing list
> >> R-sig-mixed-models@r-project.org
> >>
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> >>
> >>
> >> End of R-sig-mixed-models Digest, Vol 26, Issue 33
> >> **************************************************
> >
> >
> >
> >
> > _______________________________________________
> > R-sig-mixed-models@r-project.org mailing list
> >
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> >
> 
> _______________________________________________
> R-sig-mixed-models@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models


      


	[[alternative HTML version deleted]]


