Reputation: 980
I am trying to do chi square test on two categories of biological data. I have a data frame like this:
Brain, Cerebelum, Heart, Kidney, liver, testis
expected 3 66 1 44 34 88
observed 6 57 4 45 35 69
structure(list(Brain = c(3L, 6L), Cerebelum = c(66L, 57L), heart = c(1L,
4L), kidney = 44:45, liver = 34:35, testis = c(88L, 69L)), .Names = c("Brain",
"Cerebelum", "heart", "kidney", "liver", "testis"), class = "data.frame", row.names = c("rand",
"cns"))
I did the test using Python:
from scipy.stats import chisquare
chisquare(obs,f_exp=exp)
which gives result as:
Power_divergenceResult(statistic=17.381684491978611, pvalue=0.0038300192430189722)
I tried to replicate the results using R, so I made the csv file, imported to R as dataframe and run the code as:
d<-read.csv(file)
chisq.test(d)
Pearson's Chi-squared test
data: d
X-squared = 4.9083, df = 5, p-value = 0.4272
why the chi squared value and P value is different in python and R?, As I calculated by hand using the simple (O-E)^2/E formula, the chi square value is equal to 17.38 as calculated by python but I can not figure out how R calculate the value of 4.90.
Upvotes: 1
Views: 2419
Reputation: 2743
I can answer your first question.
chisq.test
, when you give it a matrix with > 2
rows and columns, treats it as two-dimensional contingency table and tests for independence between observations along the rows and columns. Here's an example and another one.
scipy.stats.chisq
on the other hand just does the X = sum( (O_i-E_i)^2 / E_i)
familiar from the definition of the test stat.
So how to square the circle? First, pass R
the observed values, then define the expected probabilities in argument p
. Second, you also need to stop R from doing a default continuity correction.
e <- d[1, ]
o <- d[2, ]
chisq.test(o, p = e / sum(e), correct = FALSE)
voila
Chi-squared test for given probabilities
data: o
X-squared = 17.139, df = 5, p-value = 0.004243
PS Tricky question for SO, possibly better for crossvalidated?
Note that R's default correction may be a good thing vs scipy
. Whether that is true is definitely for crossvalidated.
PPS
The help in ?chisq.test
is a litttttttle hard to parse, but I think this is all in there somewhere ;)
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.
If ‘x’ is a matrix with at least two rows and columns, it is taken
as a two-dimensional contingency table: the entries of ‘x’ must be
non-negative integers. Otherwise, ‘x’ and ‘y’ must be vectors or
factors of the same length; cases with missing values are removed,
the objects are coerced to factors, and the contingency table is
computed from these. Then Pearson's chi-squared test is performed
of the null hypothesis that the joint distribution of the cell
counts in a 2-dimensional contingency table is the product of the
row and column marginals.
and
correct: a logical indicating whether to apply continuity correction
when computing the test statistic for 2 by 2 tables: one half
is subtracted from all |O - E| differences; however, the
correction will not be bigger than the differences
themselves. No correction is done if ‘simulate.p.value =
TRUE’.
Upvotes: 2