R's chisq.test degrees-of-freedom compared to manual calculation

Question

I'm confused about the degrees-of-freedom discrepancy between R's chisq.test function, and doing the calculation "by hand".

Using the code below, the test statistics in the two approaches agree exactly but chisq.test reports 11 degrees-of-freedom, whereas doing it "by hand" I would have used k-p-1=10 where k=number of categories (12), p = number of parameters estimated from the data (1 on account of estimating the mean). I don't see how the chisq.test can know the number of parameters estimated from the data and how the number of degrees-of-freedom can agree?

k <- 0:11
obs <- c(57,203,383,525,532,408,273,139,45,27,10,6)

#Get mean...
x <- rep(k,times=obs)
lambda <- mean(x)

#probabilites under H0 (Poisson model) and expected values ek...
p <- dpois(0:10,lambda)
# p(k >= 11)
prest <- 1-sum(p)
pall <- c(p,prest)
ek <- pall*sum(obs)

#R test...
xsq <- chisq.test(obs,p=pall)

#Manual...
xsqts <- sum((obs-ek)**2/ek)
xsqdf <- length(k)-1-1
pval <- 1-pchisq(xsqts,xsqdf)
```

Answer 1

Besides the points highlighted in the comments, I guess the issue is with the way degrees of freedom is calculated in chisq.test when x (here, obs) is a vector. There is an explanation in the documentation:

If x is a matrix with one row or column, or if x is a vector and y is not given, then a goodness-of-fit test is performed (x is treated as a one-dimensional contingency table). The entries of x must be non-negative integers. In this case, the hypothesis tested is whether the population probabilities equal those in p, or are all equal if p is not given.

This wikipage explains how to calculate the degrees of freedom for a goodness-of-fit test:

For a test of goodness-of-fit, df = Cats − Parms, where Cats is the number of observation categories recognized by the model, and Parms is the number of parameters in the model adjusted to make the model best fit the observations: The number of categories reduced by the number of fitted parameters in the distribution.

So, in this case, it should be length(k)-1, 11. And this is the number reported by chisq.test.

Edit:

Snedecor and Cochran (1989, p.77) explains the calculation:

The theoretical frequency distribution of $s^2$ under normality depends on the population variance $\sigma^2$ and on the degrees of freedom $\nu$ in $s^2$. Since both of these were known when the theoretical distribution in table 5.9.1 [a table with frequencies and information about $\sigma^2$] was fitted, we say that the theoretical distribution is completely specified. In this case, the number of degrees of freedom in $\chi^2$ is $(k — 1)$, where $k$ is the number of classes whose contributions are summed in finding $\chi^2$- Note that the degrees of freedom in $\chi^2$ do not involve the sample size $n$ [emphasis mine].

In the next paragraph, they continue:

If we want to test whether this theoretical distribution fits but do not know the value of $\sigma^2$, we can use the average of the 511 sample variances $s^2$ [in the above mentioned table] as an estimate of $\sigma^2$ when fitting. The number of degrees of freedom in $\chi^2$ is then $(k —2)$…. Similarly, in fitting a normal distribution, the degrees of freedom in $\chi^2$ are $(k — 1)$ if $\mu$ and $\sigma$ are both known for the populations when fitting,…,but the degrees of freedom are $(k — 3)$ if $\bar X$ and $s$ from the sample are used as estimates of $\mu$ and $\sigma$ when fitting. The rule is

$$\text{df = (k - 1 - number of fitted parameters)}$$

So, answer to your last question (as @Glen_b pointed out): it can’t. For example, if you use pearson.test from nortest package you can see that it calculates ($k-3$) degrees of freedom (if you specify n.classes and adjust=TRUE) because sample is used to estimate parameters (but there is a warning, see the note). This is not the case for chisq.test. [I hope these long quotes don't make the answer difficult to read.]