Proper contingency table analysis for old diabetes study

I was chasing down a citation about the incidence of postoperative wound infection in diabetic patients, and I found this 1998 abstract by Pomposelli with the following statement:

"In patients with hyperglycemia (> 220 mg/dL) on POD 1, the infection rate was 2.7 times that observed (31.3% vs 11.5%) in diabetic patients with all serum glucose values < 220 mg/dL. When minor infection of the urinary tract was excluded, the relative risk for "serious" postoperative infection increased to 5.7 when any POD 1 blood glucose level was > 220 mg/dL."

I pulled the article and the contingency table data is as follows (excerpted):

Table III: Blood glucose and infection rate contingency tables
POD 1    Highest glucose    Infected     Uninfected    p value
<=220 mg/dL        3            23            .05
>220 mg/dL        21           46

Table IV: Postoperative day 1 blood glucose and infection rate
contingency table excluding urinary tract infections
Highest glucose    Infected     Uninfected    p value
<=220 mg/dL        1            23            .03
>220 mg/dL        15           46


The caption for both says "Actual p values of chi-squared shown."

When I run this in R I get:

> pod1AllInfections <- data.frame(infected=c(3,21),uninfected=c(23,46))
> row.names(pod1AllInfections) <- c("above220","below220")
> chisq.test(pod1AllInfections,simulate.p.value=TRUE)

Pearson's Chi-squared test with simulated p-value (based on 2000 replicates)

data:  pod1AllInfections
X-squared = 3.8372, df = NA, p-value = 0.06997


and for the other table:

> pod1ExcludeUTIs <- data.frame(infected=c(1,15),uninfected=c(23,46))
> row.names(pod1ExcludeUTIs) <- c("above220","below220")
> chisq.test(pod1ExcludeUTIs,simulate.p.value=TRUE)

Pearson's Chi-squared test with simulated p-value (based on 2000 replicates)

data:  pod1ExcludeUTIs
X-squared = 4.7017, df = NA, p-value = 0.03598


My question is, am I doing the chi-square calculation correctly in R? I don't understand the simulate.p.value=TRUE flag that well despite reading what I can on Google.

My second question is, for this study what's the correct contingency table analysis to do?

My last question is about Table IV -- is it biostatistically kosher to just exclude a subgroup of patients (in this case, those with UTIs) from the table?

-

There is not a right answer here. You have something which appears to potentially have a substantial impact (see the sample estimates of the odds ratios at the end) but is on the edge of a 95% statistical significance. Excluding a group of patients may be valid if it was planned before the data was taken, but doing it after seeing the data is more dubious.

It is quite easy to reproduce the original $p$ values using R

> chisq.test(matrix(c(3,21,23,46), ncol=2), correct=FALSE)

Pearson's Chi-squared test

data:  matrix(c(3, 21, 23, 46), ncol = 2)
X-squared = 3.8372, df = 1, p-value = 0.05013

> chisq.test(matrix(c(1,15,23,46), ncol=2), correct=FALSE)

Pearson's Chi-squared test

data:  matrix(c(1, 15, 23, 46), ncol = 2)
X-squared = 4.7017, df = 1, p-value = 0.03013

Warning message:
In chisq.test(matrix(c(1, 15, 23, 46), ncol = 2), correct = FALSE) :
Chi-squared approximation may be incorrect


The warning probably comes because the "expected" value in the table IV for Infected <=220 mg/dL is about $4.52$, i.e. less than $5$.

Do the continuity correction and you get less significant $p$ values

> chisq.test(matrix(c(3,21,23,46), ncol=2), correct=TRUE)

Pearson's Chi-squared test with Yates' continuity correction

data:  matrix(c(3, 21, 23, 46), ncol = 2)
X-squared = 2.8725, df = 1, p-value = 0.0901

> chisq.test(matrix(c(1,15,23,46), ncol=2), correct=TRUE)

Pearson's Chi-squared test with Yates' continuity correction

data:  matrix(c(1, 15, 23, 46), ncol = 2)
X-squared = 3.4601, df = 1, p-value = 0.06287

Warning message:
In chisq.test(matrix(c(1, 15, 23, 46), ncol = 2), correct = TRUE) :
Chi-squared approximation may be incorrect


though what this is really pointing to is a small sample size.

An alternative approach is Fisher's Exact Test, making the results more significant but not quite as much as originally

> fisher.test(matrix(c(3,21,23,46), ncol=2) )

Fisher's Exact Test for Count Data

data:  matrix(c(3, 21, 23, 46), ncol = 2)
p-value = 0.06519
alternative hypothesis: true odds ratio is not equal to 1
95 percent confidence interval:
0.05007644 1.12188507
sample estimates:
odds ratio
0.2890619

> fisher.test(matrix(c(1,15,23,46), ncol=2) )

Fisher's Exact Test for Count Data

data:  matrix(c(1, 15, 23, 46), ncol = 2)
p-value = 0.03294
alternative hypothesis: true odds ratio is not equal to 1
95 percent confidence interval:
0.003041998 0.989829479
sample estimates:
odds ratio
0.1355669

-
+1, this is a really good answer & there is a lot of info here, but 1 suggestion: you make a subtle point in your first sentence, but then move on to more concrete issues. It might be worthwhile to have a little semi-philosophical explication of what it means that "there is not a right answer here". –  gung Apr 7 '12 at 15:52
@gung: I wrote the fist paragraph last as a postscript and then moved it up. The issue of small sample sizes, publication bias and related issues are important, but probably rather wider than this particular question. –  Henry Apr 7 '12 at 16:34