By the phrase by hand
, I of course don't mean completely,
what I mean is that the computer doesn't do everything for you.
In this section, we just use R as a calculator that has built-in chi-square tables.
chisq.test
More or less self-explanatory. For more info
see the
on-line help for the chisq.test
function.
The main point that can be found out from looking at the help is that
the default null hypothesis (what is used in the example above)
is all categories equal probability. Specify a different null hypothesis
with a vector p
of probabilities
of the same length as x
.
If the sample size is too small for asymptotic approximation the function warns about this and can use simulation to get an accurate P-value. Oops! No it doesn't. Silly me, expecting the documentation to actually be correct. What it actually does is use simulation (if asked) only for chi-square tests of homogeneity and independence in two-way tables (the subject of Sections 9.3 and 9.4 in DeGroot and Schervish).
When there is a compound null hypothesis, you just do the same thing except you have to
For these data the goodness of fit test fails. We don't believe these look like the results of 10,000 rolls of true dice.
I theorize that these are six-ace flats, dice shaved down on the six and ace (one) sides so that the other sides have smaller area and consequently smaller probability. This corresponds to a probability model with one free parameter. The six-ace sides have one probability, call it θ, and the other four have another probability. In order for the other sides to add to one, we must have the other probability be (1 - 2 θ) / 4.
I claim the data in the file given by the dataset URL
below and shown in the stem-and-leaf plot are Poisson.
Since we can't have an infinite set of categories for the chi-square
test we have to lump some, say we have a 7+
category.
Oops! Only three counts in the smallest category violates
the at least 5 in each cell
rule of thumb. Let's pool
the bottom two cells too.
test we have to lump some, say we have a 7+
category.
The reason the preceding stuff is called bogus
is that it
doesn't use efficient estimation. It should actually use the likelihood
for the binned data, that is for y
not x
.
The relevant bits are the MLE based on x
Rweb:> lambda.hat <- mean(x) Rweb:> print(lambda.hat) [1] 3.53 Rweb:> p.hat <- fred(lambda.hat) Rweb:> print(p.hat) 1- 2 3 4 5 6 7+ 0.13275127 0.18258281 0.21483911 0.18959551 0.13385443 0.07875102 0.06762584
and the MLE based on y
Rweb:> lambda.twiddle <- out$estimate Rweb:> print(lambda.twiddle) [1] 3.513539 Rweb:> p.twiddle <- fred(lambda.twiddle) Rweb:> print(p.twiddle) 1- 2 3 4 5 6 7+ 0.13446416 0.18388611 0.21536369 0.18917219 0.13293278 0.07784408 0.06633699
almost no difference. It might make a difference if the P-value
was near the borderline of statistical significance. In this example,
where P = 0.69, it is just pedantic silliness to worry about
whether we use lambda.hat
or lambda.twiddle
.
(We only did it to show how.)
Let's try Example 9.3.1 in DeGroot and Schervish, first by hand,
then using chisq.test
. Hmmmmm. No, changed my mind.
Other way round.
chisq.test
We can't use the dataset URL
mechanism in Rweb for two-way tables,
because it won't read a table properly, just a data frame. So we enter
the data the hard way. The function rbind
binds rows together
to make a matrix.
Of course, the dimnames(x)
assignment is completely
unnecessary. It only affects the appearance of the print(x)
statement. The calculation only depends on the numbers.
No real reason for this, just showing that the chisq.test
function does do the same thing as calculating by hand.
If you only want to see the expected counts, the chisq.test
function does that. No need for hand calculation.
Consider the following data, a 3 by 3 table
0 9 8 8 8 11 7 19 30
Under the multiplicative model, there will not be at least 5 expected in
each cell and the chisq.test
function warns about this.
However we can get an accurate P-value by simulation
There's no need to use a simulation size as small as the default. The person who wrote this function must have had a really slow computer.