4
$\begingroup$

I have got an answer for it from Spacedman. But I am not entirely satisfied with the answer as it does not give me any sort of value (p or z value). So I am re-framing my question and posting it again. No offences to Mr.Spacedman.

I have a dictionary of say 61000 elements and out of this dictionary, I have two sets. Set A contains 23000 elements and Set B contains 15000 elements and an overlap of Set A and Set B gives 10000 elements. How can I estimate a p-value or z value to show that this overlap is significant and is not occuring by chance or vice versa.

What I have been suggested till now includes MonteCarlo simulation methods. Is it possible to have an analytical method.

Thank you in advance.

Improve this question
$\endgroup$
7
  • $\begingroup$ Significant of what? Suppose A is {a,b,c,d,e,f} (6 elements) and B is {d,e,f,g,h,i,j} (7 elements). The overlap of A and B is {d,e,f} (3 elements). So what? What are you asking? $\endgroup$ – Spacedman Aug 22 '11 at 13:27
  • $\begingroup$ these elements are genomic co-ordinates, and thus overlap may have been by chance. Hence some sort of p-value either by phyper (from R) or Fischer test is needed. $\endgroup$ – Angelo Aug 22 '11 at 13:35
  • $\begingroup$ So you're looking at a partition of (23706+14557-10752) unique elements into two sets of size 23706 and 14557 and want the distribution of the size of A intersection B assuming random sampling without replacement? $\endgroup$ – Spacedman Aug 22 '11 at 14:05
  • $\begingroup$ Why do you need a p or z value? $\endgroup$ – richiemorrisroe Aug 23 '11 at 8:36
  • 1
    $\begingroup$ To mods or anyone more familiar with CV: surely this site has a FAQ for contingency tables - if so, then it would be best to link to that, and perhaps close as a duplicate. So far, I cannot seem to find a good, clean example question here regarding contingency tables. $\endgroup$ – Iterator Aug 24 '11 at 20:42
8
$\begingroup$

If my comment is right, then you can Monte-carlo simulate it:

sim=unlist(lapply(1:10000,
 function(i){A=sample(1:27511,23706);B=sample(1:27511,14557);return(sum(A %in% B))}))
hist(sim)

Probably neater ways to do that loop but whatever.

Your 10752 is waaaay over to the left of my histogram, so significantly fewer common elements than expected by chance.

There may be some exact test that does the same. There's probably a Normal approximation - in which case it looks about 20 sigma off the mean:

hist(sim,xlim=c(10752,12660))
Improve this answer
$\endgroup$
8
  • $\begingroup$ Thank you, for your reply. Let elaborate a bit first. I have a dictionary from which two lists are drawn and I have to compute the significance of overlap of these two lists whether it is is significant or is it happening by chance. With this given criteria can I still use your command list for MC? looking for a reply. $\endgroup$ – Angelo Aug 22 '11 at 14:18
  • 1
    $\begingroup$ Pretty much. My 'dictionary' is the vector of numbers 1:27511. I take two samples without replacement - one of size 23706 and one of size 14557. I then compute the size of the overlap. The 27511 is the 23706+14557 (A+B) minus 10752 (A&B) because otherwise you count the intersection twice. Draw a Venn diagram. $\endgroup$ – Spacedman Aug 22 '11 at 14:26
  • $\begingroup$ Thank you very much it works very nicely, could please also think of a way to get a value which (p value or z value) could be an indicative of the overlap being significant or not. $\endgroup$ – Angelo Aug 22 '11 at 14:48
  • $\begingroup$ Standard Monte-Carlo test stuff - if you do 99 simulations and your data is still the largest (or smallest) then p<0.01. Until you've got simulations beyond your data value you can only quote a 'p less-than' value. Unless you maybe do some smoothing or something kludgy... $\endgroup$ – Spacedman Aug 22 '11 at 15:10
  • $\begingroup$ Thank you again. Is their a way to be in touch with you other than this forum (not to bother you with statistical question every now and then but for expert advice and comments). Further, I just noticed I cannot accept your answer as I do not have 15 points. :( But your answer is excellent. $\endgroup$ – Angelo Aug 22 '11 at 15:18
7
$\begingroup$

A model for this situation is to put 61000 ($n$) balls into an urn, of which 23000 ($n_1$) are labeled "A". 15000 ($k$) of these are drawn randomly without replacement. Of these, $m$ are found to be labeled "A". What is the chance that $m \ge 10000$?

The total number of possible samples equals the number of $k$-element subsets of an $n$-set, $\binom{n}{k}$. All are equally likely to be drawn, by hypothesis. Let $i \ge 10000$. The number of possible samples with $i$ A's is the number of subsets of an $n_1$-set having $i$ A's, times the number of subsets of an $n-n_1$-set having $k-i$ non-A's; that is, $\binom{n_1}{i}\binom{n-n_1}{k-i}$. Summing over all possible $i$ and dividing by the chance of each sample gives the probability of observing an overlap of $m = 10000$ or greater:

$$\Pr(\text{overlap} \ge m) = \frac{1}{\binom{n}{k}} \sum_{i=m}^{\min(n_1,k)} \binom{n_1}{i}\binom{n-n_1}{k-i}.$$

This answer is exact. For rapid calculation it can be expressed (in closed form) in terms of generalized hypergeometric functions; the details of this expression can be provided by a symbolic algebra program like Mathematica. The answer in this particular instance is $3.8057078557887\ldots \times 10^{-1515}$.

We can also use a Normal approximation. Coding A's as 1 and non-A's as 0, as usual, the mean of the urn is $p = 23000/61000 \sim 0.377$. The standard deviation of the urn is $\sigma = \sqrt{p(1-p)}$. Therefore the standard error of the observed proportion, $u = 10000/15000 \sim 0.667$, is

$$se(u) = \sigma \sqrt{(1 - \frac{15000-1}{61000-1})/15000} \sim 0.003436.$$

(see http://www.ma.utexas.edu/users/parker/sampling/woreplshort.htm). Thus the observed proportion is $z = \frac{u - p}{se(u)} \sim 84.28$ standard errors larger than expected. Obviously the corresponding p-value is low (it computes to $1.719\ldots \times 10^{-1545}$). Although the Normal approximation is no longer very accurate at such extreme z values (it's off by 30 orders of magnitude!), it still gives excellent guidance.

Improve this answer
$\endgroup$
5
  • $\begingroup$ No offense, but this seems to be a peculiar approach to a contingency table. $\endgroup$ – Iterator Aug 24 '11 at 20:34
  • $\begingroup$ @Iterator Is any part of it inaccurate or incorrect? How does the approximate calculation compare to the calculations involved in a chi-squared test? :-) Is there any ambiguity about the assumptions, the probability model, or what is being estimated? $\endgroup$ – whuber♦ Aug 24 '11 at 21:12
  • $\begingroup$ It is a correct solution but I think it is for the wrong problem. The difference wrt to the contingency table is that the null hypothesis of the contingency table is well defined relative to the independent probabilities of the bins. This approach tests a different statistic; relating it to null hypotheses will be somewhat more cumbersome. However, a it is certainly feasible to describe a contingency table of counts in an urn setting. $\endgroup$ – Iterator Aug 24 '11 at 21:47
  • $\begingroup$ @Iterator I'm afraid I don't follow. What exactly is the "wrong problem"? Are you saying you disagree with this model of the situation? What do you think is the right model? And what's "cumbersome" about the null hypothesis? (Please recall that the chi-squared test is purely an approximation that is justified solely by an appeal to a model like the one I started with.) $\endgroup$ – whuber♦ Aug 25 '11 at 3:19
  • $\begingroup$ No problem: suppose you set up the problem as 4 urns along 2 axes, and tossed balls in. The total could be 61K balls, and the question can be rephrased as to whether or not the axes are independent. In contrast, your model imposes a limit on axis 2 (the second category, i.e. the smaller # of balls) as a subset of the # chosen for the larger set. This peculiar and wouldn't generalize (i.e. why does one have to choose the larger set?). As for the right answer, the chi-square p-value is also infinitesimal. $\endgroup$ – Iterator Aug 25 '11 at 3:57
6
$\begingroup$

You have a large enough sample size that a Chi-square test is reasonable. If you are using R, then chisq.test() is the function to utilize; its help page can be found via ?chisq.test.

Improve this answer
$\endgroup$
3
  • $\begingroup$ would you mind sparing a minute and writing a sample syntax for my dataset. It would be of great help. $\endgroup$ – Angelo Aug 23 '11 at 10:12
  • $\begingroup$ The help for chisq.test will show you how to use it. Even so, you will find R saying things like "p-value < 2.2e-16" with data way off the mean. $\endgroup$ – Spacedman Aug 23 '11 at 11:34
  • $\begingroup$ @Angelo: You should try it out. Don't be afraid to look at the help and plug your own data in. It won't bite. :) You should also know what a chi-squared test actually is, rather than just select an arbitrary test someone told you about on the internet. $\endgroup$ – Iterator Aug 23 '11 at 14:34
5
$\begingroup$

Well, if nothing else, you can do it by simulation:

Draw 23000 random elements out of 61000 and draw 15000 random elements from the same 61000. Now count the number of overlapped items.

Repeat 100000 times (should not take all too long): now you have an empirical distribution of the number of overlapped items, and you can easily find an empirical p-value for 10000 elements.

Code like the following will do this:

countOverlap<-function(total=61000, numgA=23000, numgB=15000, replace=FALSE){
    groupA<-sample.int(total, numgA, replace=replace)
    groupB<-sample.int(total, numgB, replace=replace)
    return(length(intersect(groupA, groupB)))
}

tmpres<-replicate(1000, countOverlap(total=100, numgA=20, numgB=30))

#if your true observed value was 10:
pval<-mean(tmpres >= 10)
pval

However, it is possible that an analytical solution exists, and that simulation is not needed. Maybe someone else can provide that for you.

Improve this answer
$\endgroup$
6
  • $\begingroup$ how about changing the last part (where you calculate pval)to something like this: pnorm(10,mean=mean(tmpres), sd=sd(tmpres), lower.tail=T). Please give your comment. $\endgroup$ – Angelo Aug 24 '11 at 14:50
  • $\begingroup$ If you have some reason to assume these numbers of overlaps to be normally distributed, it can be used. However, especially since they are discrete numbers (and normal is a continuous distribution), this seems like a bad idea to me (and, if your number of repeats in the above is high enough, should not return more realistic p-values) $\endgroup$ – Nick Sabbe Aug 24 '11 at 15:16
  • $\begingroup$ Yes, I checked these numbers and they are normally distributed. Mean=Median=Mode (they are almost equal). Well the problem with your pval calculation is, either I am getting 1 or 0. Hence, what will be a good trade off? $\endgroup$ – Angelo Aug 24 '11 at 15:24
  • $\begingroup$ Having Mean=Median=Mode is not the same as having a normal distribution! I expect you to get very small p-values, because your situation is extreme (there are way more overlapped items than expected by chance - as also indicated by @whuber). Typically this is reported (in research papers) as "p-value << 0.001" or similar. $\endgroup$ – Nick Sabbe Aug 24 '11 at 15:29
  • $\begingroup$ Yes I agree with you. Furthermore, I am gettingvery small p-values. Hence, the (best) approach is to plot a histogram and check whether my observed value is towards left or right of histogram (as suggested by Spacedman). Do you agree to this approach? $\endgroup$ – Angelo Aug 24 '11 at 15:39
0
$\begingroup$

The correct approach depends on what is fixed. If the size of the 2 sets are "fixed" then this is just done by a chi-squared test of independence where the 2x2 table is "in /not in" each set. If only the total number of items is fixed, then the chi-squared test is not correct.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.