this is my first asking question here... forgive me if my question is not clear enough.

I have two datasets; one is ground truth dataset and another is perturbed version of the same dataset. I want to test if a regression model (algorithm) behaves similarly on two datasets. More specifically, I want do things like following:


fit1 <- mgcv::gam(y ~ s(x), family = poisson, data = my_data)
fit2 <- mgcv::gam(y ~ s(x), family = poisson, data = my_data.perturbed)
# should return p-value for if two models are same
test_two_model_are_similar(fit1, fit2)

This is analogous to comparing coefficients of two parametric models (as in here and here). In parametric models, "test if two models behave in similar way" can be defined as test that coefficients of the models are statistically indistinguishable (I do not need to statistically prove that two models are indistinguishable; I only want to show that I cannot reject the null hypothesis that two models are different). However, I want to do similar things on a non-parametric model such as GAM with splines, in which case I cannot compare coefficients (or can I?).

I do not intend to compare the goodness of fit; so doing anova(fit1, fit2) or comparing BIC is not appropriate.

So, my question can be broken down into following questions:

  • Is there any direct way (preferably with R) that I can do the above stuffs?
  • If I can sample $y_i$ given $x_i$ from model 0 (so I get $P_0(y_i|x_i)$) and calculate probability of $P_1(y_i|x_i)$, I think I can numerically estimate the overlaps between outputs of the two models. Is there any way on R that:
    • I can (efficiently) sample $x_i$ to feed into model 0?, and
    • calculate probability of $y_i$ given $x_i$ given model 1?
  • $\begingroup$ One option might be to model s(x) for both data sets in a single model wherein we model s(x) for data as the reference level and then s(x) for data.perturbed` as a smooth difference between the original s(x) and the perturbed data. This is analogous to how R treats factors with the reference level being the intercept and the other coefs for the other levels represent differences between each level and the reference level. Then we have a Wald-like test for the second (difference) smooth that is analogous to testing s(x)data.perturbed == 0. Does that sound useful? $\endgroup$ – Gavin Simpson Mar 5 '19 at 4:09
  • 4
    $\begingroup$ @koreyou Your problem is that it is likely that 'data set A' and 'perturbed version of data set A' cannot be treated either as the same data nor as independent data, but it is not clear from the above how to treat their dependence or even whether the second can be treated as a random variable. Note further that "the two models are different" is not a suitable null hypothesis. You will need to explain more about what you're doing. $\endgroup$ – Glen_b Mar 5 '19 at 9:34