1
$\begingroup$

I have a data set, whose elements (say 1,000 elements) are numeric with values between 0 and 1. It is safe to think the elements are independent from each other, i.e, the order of these elements does not matter at all. Also, from the histogram of the data set, I can tell it does not follow any parametric distribution (normal or others).

I used a modeling tool (not in R) to simulate 100 realizations of this data set given a set of parameters. In each realization, the number of elements is not necessarily equal to 1,000.

My goal is simple: Does the model fit the data well?

I plotted the histograms of the simulated data against the empirical data one (using both counts and probability density) and was sort of able to tell the quality of the fit, but would like to have a solid statistical way to quantitatively reach my goal (say showing the p-value of this fit is > 10% or so).

Is there a good way in R or Python I can do so? I looked up the standard R fitting methods, such as lm() and glm() in R, but they are for fitting a line to the data, not doing something I would like to have. I believe that this is simply because I am not familiar with R and statistics given this question does not seem to be weird and difficult.

Thanks, Ben

| cite | improve this question | |
$\endgroup$
  • 4
    $\begingroup$ ?ks.test in R (Kolmogorov-Smirnov) $\endgroup$ – Ben Bolker Feb 24 '12 at 20:39
  • $\begingroup$ Your question seems to contradict itself because you say your original data does not follow a parametric distribution, but then you go on to say somehow you have simulated similar data using a set of parameters. Probably what you mean is it does not follow any simple parametric distribution, but perhaps it can be seen as a mixture eg of two or three different beta distributions? $\endgroup$ – Peter Ellis Feb 25 '12 at 1:09
  • $\begingroup$ @Peter.. My so-called simulation data were obtained from a partial differentiatial equation based modeling approach based on some a priori knowledge (parameters). My goal is not to estimate these parameters, but to tell if the model fits the data well (quantitatively). I have no idea what the true distribution of the data (either empirical or simulation) would be. $\endgroup$ – Ben Feb 25 '12 at 3:16
6
$\begingroup$

You can use the qqplot function like @Peter Flom suggests as one way to compare visually. Another visual approach is to create a set of graphs where one is the real data and the others are simulated data, but you don't know which is the real data (or a sample of the real data), then see if you can pick out the real data. The vis.test function in the TeachingDemos package for R helps with this. This gives more of a close enough test than an exact match test.

The ks.test function in R can take 2 vectors of numbers and tests the null hypothesis that they come from the same distribution. Note however that with large sample sizes this test can find significance for differences so small that most would not care about and for small sample sizes it may not have the power to differentiate between distributions that we would care about the difference of.

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ @Greg... Thanks! I will take a look at the ks.test and vis.test. $\endgroup$ – Ben Feb 24 '12 at 21:56
5
$\begingroup$

Let me see if I understand what you've got: One set of "real data" with 1000 values and 100 sets of "model data" each with about 1000 values.

You could combine all 100 sets and then, in R, run something like qqplot(real, fake).

Or you could do 100 such plots (one for each set).

I don't think the concept of a p-value makes much sense here, as you don't have a null hypothesis that you are attempting to reject, or a test-statistic that you expect not to be large given a null.

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ @Peter...I think my null hypothesis is that the real and model data come from the same (though unknown) distribution. A p-value (or any statistical evidence) can give a solid quantitative judgement on the fit, rather than simply rely on a graphically qualitative judgement. :) $\endgroup$ – Ben Feb 24 '12 at 21:46
  • $\begingroup$ @Peter... most importantly, thanks for the tips on qqplot. surely will give it a try! $\endgroup$ – Ben Feb 24 '12 at 21:55
  • 2
    $\begingroup$ Unfortunately, Ben, when you create the model in response to your review of the data, then--unless you are really bad at modeling--any p-values for goodness of fit you subsequently compute will be hugely over-optimistic. One way around this problem is to retain some of the original data for confirmatory testing later, but in this case that horse left the barn long ago... $\endgroup$ – whuber♦ Feb 27 '12 at 16:46
  • $\begingroup$ Ben if you want something "solid" or "quantitative" then define "well" in your original question in a solid, quantitative way. Even beyond the valid point that @whuber raised, there is no case when a p-value tells you if a model fits "well". A p-value tells you if a test statistic as extreme as that you get in a sample is likely if the null hypothesis is true in the population. There are very few cases where this question is of much interest, and this isn't one of them, even if you could define it. $\endgroup$ – Peter Flom♦ Feb 27 '12 at 19:20
0
$\begingroup$

For example if you have data x<-c(1,2,3,5,4,6,9) and you have so many distributions to see which one is best. you can do it in R for example, I will consider first Poisson distribution.

dburr=function(obs, c =2.1000001) exp(-c)*c^x/factorial(x)
library(FAdist)
set.seed(45)
#obs <- c(1,1,0,0,0,8,0,5,0,0)
#obs<-c(rep(0,17),rep(1,9),rep(2,7),rep(3,3),4,5,6,6)
obs<-c(rep(0,64),rep(1,17),rep(2,10),rep(3,9))
#View(obs)
library(MASS)
jk=fitdistr(x = obs,densfun = dburr,start = list(c=0.001), # need to provide named list of starting values
                                 lower = list(c=0.02)) # and named list of lower bounds since c, k > 0
jk$loglik
    jk$vcov
jk$estimate
    jk$sd
jk$n
AIC(jk)
BIC(jk)

Similarly define other distribution density functions in the above code and calculate AIC, BIC, log likelihood....the distribution with minimum AIC, BIC will be good as compared to one having higher value for AIC, BIC.

| cite | 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.