I have 200000 simulations and I want to use Approximate Bayesian Computation (ABC) to determine the best 1000, based on specific targets. These simulations have 12 parameters (my priors, dependent variables). For each simulation I extracted summary statistics and I have four targets.

All my parameters are on different scales, for example here is a summary of the first 5:

> summary(sample_parameters)[, 1:5]
 standard_deviation growth_rate_max  density_max     dispersal_p     dispersal_r   
 Min.   :0.000      Min.   :0.560   Min.   : 1000   Min.   :0.050   Min.   :100.0  
 1st Qu.:0.090      1st Qu.:0.940   1st Qu.: 4625   1st Qu.:0.110   1st Qu.:244.0  
 Median :0.175      Median :1.320   Median : 8250   Median :0.175   Median :387.5  
 Mean   :0.175      Mean   :1.316   Mean   : 8250   Mean   :0.175   Mean   :387.5  
 3rd Qu.:0.260      3rd Qu.:1.690   3rd Qu.:11875   3rd Qu.:0.240   3rd Qu.:531.0  
 Max.   :0.350      Max.   :2.070   Max.   :15500   Max.   :0.300   Max.   :675.0 

I am using the 'abc' package in R. I am moving beyond the rejection method, and I use a regression adjustment to account for the imperfect match between simulations and observations. Because I have many parameters, this imperfect match can be more pronounced and credible intervals from the regression approach can be inflated. I therefore use a logit transformation with logit bounds, to restrict the adjustment of the posteriors within my parameters' ranges.

> head(targets)
   fossil_occupancy_1500AD colonization_penalty_flat RMSE_flat abundance_trend_tot
1:                     720                         0         0         108.0695

> head(abc_df)
   fossil_occupancy_1500AD colonization_penalty_flat RMSE_flat abundance_trend_tot
1:                       1                     -1272 1239.9636              -8.107
2:                     353                     -1272 1283.0134            -217.652
3:                     717                      -689 1285.4261            3753.685
4:                      56                     -1272  933.9658            -256.001
5:                     278                     -1272  664.7752              84.517
6:                     289                     -1272  688.5428            -480.711

nsim <- 1000

set.seed(1234); abc_ridge <- abc(target = targets,
                                 param = sample_parameters,
                                 sumstat = abc_df,
                                 tol = nsim/200000,
                                 method = "ridge",
                                 transf = 'logit',
                                 logit.bounds = ranges)

Because all my parameters are on different scales, my first thought was to centre-scale my data before running the ABC. However, while I was reading about ABC approaches, I did not find mentions of scaling the data before using this method.

I decided to run an ABC with scaled data (logit(scale(x))) and another with the priors in their own scales (logit(x); as shown in the code above). When I look at the resulting qq-plots for the residuals, the logit(x) shows a much straighter line compared to the logit(scale(x)). However, I get a much closer and balanced match to my targets if I scale the data first.

So my question is, why centre-scaling the input data gives such different results compared to using the priors in their own scales? And is it theoretically wrong to centre-scale the input data before running an ABC analysis?

Thank you!

  • $\begingroup$ I am afraid I do not understand the question but the scaling in ABC should bear upon the summary statistics, not the raw data (unless used as summary statistics). $\endgroup$
    – Xi'an
    Sep 6 at 8:20

Your Answer

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

Browse other questions tagged or ask your own question.