Re: [R] chisq.test(): standardized (adjusted) Pearson residuals
On Aug 19, 2011, at 1:28 PM, Stephen Davies wrote: I'm using chisq.test() on a matrix of categorical data, and I see that the "residuals" attribute of the returned object will give me the Pearson residuals. That's cool. However, what I'd really like is the standardized (adjusted) Pearson residuals, which have a N(0,1) distribution. Is there a way to do that in R (other than by me programming it myself?) ?scale -- David Winsemius, MD West Hartford, CT __ R-help@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.
Re: [R] chisq.test(): standardized (adjusted) Pearson residuals
On Aug 19, 2011, at 20:40 , David Winsemius wrote: > > On Aug 19, 2011, at 1:28 PM, Stephen Davies wrote: > >> I'm using chisq.test() on a matrix of categorical data, and I see that the >> "residuals" attribute of the returned object will give me the Pearson >> residuals. >> That's cool. However, what I'd really like is the standardized (adjusted) >> Pearson residuals, which have a N(0,1) distribution. Is there a way to do >> that >> in R (other than by me programming it myself?) > > ?scale chisq.test(...)$stdres, more likely. -- Peter Dalgaard, Professor, Center for Statistics, Copenhagen Business School Solbjerg Plads 3, 2000 Frederiksberg, Denmark Phone: (+45)38153501 Email: pd@cbs.dk Priv: pda...@gmail.com "Døden skal tape!" --- Nordahl Grieg __ R-help@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.
Re: [R] chisq.test(): standardized (adjusted) Pearson residuals
On Aug 20, 2011, at 3:43 AM, peter dalgaard wrote: On Aug 19, 2011, at 20:40 , David Winsemius wrote: On Aug 19, 2011, at 1:28 PM, Stephen Davies wrote: I'm using chisq.test() on a matrix of categorical data, and I see that the "residuals" attribute of the returned object will give me the Pearson residuals. Actually they are not an attribute in the R sense, but rather a list value. That's cool. However, what I'd really like is the standardized (adjusted) Pearson residuals, which have a N(0,1) distribution. Is there a way to do that in R (other than by me programming it myself?) ?scale chisq.test(...)$stdres, more likely. Agree that does have a much greater chance of keeping the questioner in the mainstream of statistics terminology and is most likely what he was looking for, but do not think the result will in general have an N(1,0) distribution. I believe the correct statement is that standardized residuals would (in the statistical "asymptotic" sense) have an N(1,0) distribution if and when the null hypothesis of marginal homogeneity were true, but should not be N(1,0) in any case when an alternate hypothesis holds. My error was in taking the questioner's request at face value. -- David Winsemius, MD West Hartford, CT __ R-help@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.
Re: [R] chisq.test(): standardized (adjusted) Pearson residuals
> >>> I'm using chisq.test() on a matrix of categorical data, and I see > >>> that the > >>> "residuals" attribute of the returned object will give me the > >>> Pearson residuals. > > Actually they are not an attribute in the R sense, but rather a list > value. Oh. I was just going by: > attributes(my.chisq.test) $names [1] "statistic" "parameter" "p.value" "method""data.name" "observed" [7] "expected" "residuals" $class [1] "htest" which I interpreted as "this object has 8 attributes, called 'statistic', 'parameter', ..., 'residuals'." Is that not the right terminology? > >>> That's cool. However, what I'd really like is the standardized > >>> (adjusted) > >>> Pearson residuals, which have a N(0,1) distribution. Is there a > >>> way to do that > >>> in R (other than by me programming it myself?) > >> > >> ?scale > > > > chisq.test(...)$stdres, more likely. "scale" is not what I want. As for "$stdres," that would be wonderful, but as you can see from the above list of attributes, it's not one of the 8 returned. What am I missing? - Stephen __ R-help@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.
Re: [R] chisq.test(): standardized (adjusted) Pearson residuals
On Aug 20, 2011, at 12:04 PM, Stephen Davies wrote: I'm using chisq.test() on a matrix of categorical data, and I see that the "residuals" attribute of the returned object will give me the Pearson residuals. Actually they are not an attribute in the R sense, but rather a list value. Oh. I was just going by: attributes(my.chisq.test) $names [1] "statistic" "parameter" "p.value" "method""data.name" "observed" [7] "expected" "residuals" $class [1] "htest" which I interpreted as "this object has 8 attributes, called 'statistic', 'parameter', ..., 'residuals'." Is that not the right terminology? The names attribute let's you know what characters to use if you want to access values in a list. Unless you are doing programming attributes is not a particular useful function. It is much more common to access the names attribute with the `names` function: > names(Xsq) [1] "statistic" "parameter" "p.value" "method""data.name" "observed" "expected" [8] "residuals" "stdres" So "stdres" is not an attribute but rather one value in a particular attribute called "names". You would get (much) more information by using str on the htest object as below: > str(Xsq) List of 9 $ statistic: Named num 30.1 ..- attr(*, "names")= chr "X-squared" $ parameter: Named num 2 ..- attr(*, "names")= chr "df" $ p.value : num 2.95e-07 $ method : chr "Pearson's Chi-squared test" $ data.name: chr "M" $ observed : table [1:2, 1:3] 762 484 327 239 468 477 ..- attr(*, "dimnames")=List of 2 .. ..$ gender: chr [1:2] "M" "F" .. ..$ party : chr [1:3] "Democrat" "Independent" "Republican" $ expected : num [1:2, 1:3] 704 542 320 246 534 ... ..- attr(*, "dimnames")=List of 2 .. ..$ gender: chr [1:2] "M" "F" .. ..$ party : chr [1:3] "Democrat" "Independent" "Republican" $ residuals: table [1:2, 1:3] 2.199 -2.505 0.411 -0.469 -2.843 ... ..- attr(*, "dimnames")=List of 2 .. ..$ gender: chr [1:2] "M" "F" .. ..$ party : chr [1:3] "Democrat" "Independent" "Republican" $ stdres : table [1:2, 1:3] 4.502 -4.502 0.699 -0.699 -5.316 ... ..- attr(*, "dimnames")=List of 2 .. ..$ gender: chr [1:2] "M" "F" .. ..$ party : chr [1:3] "Democrat" "Independent" "Republican" - attr(*, "class")= chr "htest" Now you can see that the values in the stdres object are really a list element and are in a table with particular row and column names. You get that object one of two ways. you ca use the "$" method as Dalgaard suggested or you can use "[[" with the name of the object: Xsq[["stdres"]] That's cool. However, what I'd really like is the standardized (adjusted) Pearson residuals, which have a N(0,1) distribution. Is there a way to do that in R (other than by me programming it myself?) ?scale chisq.test(...)$stdres, more likely. "scale" is not what I want. As for "$stdres," that would be wonderful, but as you can see from the above list of attributes, it's not one of the 8 returned. What am I missing? David Winsemius, MD West Hartford, CT __ R-help@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.
Re: [R] chisq.test(): standardized (adjusted) Pearson residuals
On Aug 20, 2011, at 18:04 , Stephen Davies wrote: > As for "$stdres," that would be wonderful, but > as you can see from the above list of attributes, it's not one of the 8 > returned. What am I missing? An upgrade, most likely. -- Peter Dalgaard, Professor, Center for Statistics, Copenhagen Business School Solbjerg Plads 3, 2000 Frederiksberg, Denmark Phone: (+45)38153501 Email: pd@cbs.dk Priv: pda...@gmail.com "Døden skal tape!" --- Nordahl Grieg __ R-help@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.
Re: [R] chisq.test(): standardized (adjusted) Pearson residuals
On Aug 20, 2011, at 12:57 PM, peter dalgaard wrote:
On Aug 20, 2011, at 18:04 , Stephen Davies wrote:
As for "$stdres," that would be wonderful, but
as you can see from the above list of attributes, it's not one of
the 8
returned. What am I missing?
An upgrade, most likely.
Whoosh. Sometimes I am simply clueless. I didn't notice that 'stdres'
was missing from the names in Stephen's output. Laura Thompson has a
very nice R/S accompaniment to Agresti's "Categorical Data Analysis"
text and she shows how to adjust the Pearson residuals to make them
"standardized". What follows is directly from pages 37-38 of her work:
#--#
resid.pear <- residuals(fit.glm, type = "pearson")
Note that the sum of the squared Pearson residuals equals the Pearson
chi-squared statistic:
sum(resid.pear^2)
[1] 69.11429
To get the standardized residuals, just modify resid.pear according to
the formula on p. 81 of Agresti.
ni<-rowSums(table.3.2.array) # row sums
nj<-colSums(table.3.2.array) # column sums
n<-sum(table.3.2.array) # total sample size
resid.pear.mat<-matrix(resid.pear, nc=3, byrow=T,
dimnames=list(c("
"Bachelor or Grad"),c("Fund", "Mod", "Lib")))
n*resid.pear.mat/sqrt(outer(n-ni,n-nj,"*") ) # standardized Pearson
residuals
FundMod Lib
You can also look at the code (once you upgrade) and the method in R
is quite similar, although the R codes calcualtes the stdres values
separately rather than adjusting the Pearson residuals
--
Peter Dalgaard, Professor,
Center for Statistics, Copenhagen Business School
David Winsemius, MD
West Hartford, CT
__
R-help@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.
