R : Calculate a P-value of a random distribution [duplicate]

This question already has an answer here:

I want to get the P-value of two randomly distributed observations x and y, for example :

> set.seed(0)
> x <- rnorm(1000, 3, 2)
> y <- rnorm(2000, 4, 3)

or:

> set.seed(0)
> x <- rexp(50, 10)
> y <- rexp(100, 11)

let's say that T is my test-statistic defined as mean(x) - mean(y) = 0 (this is H0), the P-value is then defined as : p-value = P[T>T_observed | H0 holds]. I tried doing this :

> z <- c(x,y) # if H0 holds then x and y are distributed with the same distribution
> f <- function(x) ecdf(z) # this will get the distribution of z (x and y)

then to calculate the p-value i tried this:

> T <- replicate(10000, mean(sample(z,1000,TRUE))-mean(sample(z,2000,TRUE)))
# this is  supposed to get the null distribution of mean(x) - mean(y)
> f(quantile(T,0.05)) # calculating the p-value for a significance of 5%

obviously this doesn't seem to work, what am i missing ?

marked as duplicate by whuber♦Jul 16 '17 at 16:02

• The problem is that this test statistic is useless for the purpose and you cannot define a p-value for it unless you make strong assumptions about what the underlying distribution might mean. (In other words, you haven't adequately described $H_0$.) That's because the distribution of this statistic depends not only on the true difference in means, but also on the variances. Since you don't know the variances, you're stuck. Go back to the textbooks and read about the Student $t$ statistic: it solves this problem. – whuber Jul 15 '17 at 16:32
• I have described H0 ! H0 : the two observations have the same distribution. H0 holds when mean(x) = mean(y), and therefore, H1: the two observations x and y doesn't have the same distribution and mean(x) does not equal mean(y) – Manuel Jul 15 '17 at 16:37
• Unfortunately, your test statistic is worthless for evaluating that $H_0$, because $H_0$ is not sufficiently specific. That leaves you two choices: either modify $H_0$ or change your test statistic. – whuber Jul 15 '17 at 16:38
• well this is and old programming contest, and it states literaly this : We take an observation x and an observation y, we want to know if x and y are distributed with the same distribution, so we consider $H_0$: x and y are distributed with the same distribution and $H_1$: x and y have different distributions. We consider the statistic test mean(x) - mean(y) = 0 means that $H_0$ holds. We want to determine the p-value using the distribution function of mean(x) - mean(y), that means generating random values with Monte-carlo and then using ecdf to generate the distribution function – Manuel Jul 15 '17 at 16:46
• Fair enough. The correct solution, then, is that the p-value is always $1$! Moreover, that's a wonderfully scalable solution, because it can be computed in constant time, and there are exceptionally brief programs that implement it, such as (1) in R. – whuber Jul 15 '17 at 16:51

c(x,y) is just the entries in x followed by the entries in y. You want t.test(x, y), which returns

Welch Two Sample t-test

data:  x and y
t = -11.934, df = 2759.7, p-value < 2.2e-16
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-1.2768896 -0.9164975
sample estimates:
mean of x mean of y
2.968341  4.065034

Enter ?t.test at the R command prompt to view the documentation for this function (including how to do things like change your desired value of $\alpha$ or use a one-sided rather than two-sided test).

Note that this is using the t-test for the difference in means between X and Y. Since your other example involves a pair of exponential distributions, you might be asking for an equality of distribution test. A standard one is Kolmogorov–Smirnov:

> ks.test(x, y)

Two-sample Kolmogorov-Smirnov test

data:  x and y
D = 0.2185, p-value < 2.2e-16
alternative hypothesis: two-sided

See paragraph 4 of the Wikipedia page on Kolmogorov-Smirnov for references to some other common tests of equality of distribution.

• i know about those built in functions, thing is i want to write it by myself .. and be able to give two observations x and y and find a p-value of statistic test – Manuel Jul 15 '17 at 16:16