3.3 Minimum Dose Model

In the absence of a well-bleached sample of the same mineral derived from the same source and of equivalent age to derive a reasonable estimate of \(\sigma_b\) from (Galbraith & Roberts, 2012), one my resort to a purily statistical reasoning for a chosen \(\sigma_b\) value. The calc_MinDose() (and likewise, the calc_MaxDose() and calc_FiniteMixture() models) provide an estimate of the maximum likelihood as well as the Bayesian Information Criterion (BIC), which can be used to choose between several results when applying a range of \(\sigma_b\) values.

The following script first calculates the equivalent doses from the provided BIN file and then uses the lapply() function to iteratively apply the calc_MinDose() while progressively increasing the \(\sigma_b\) value. We then choose the model with the lowest BIC score as our final estimate.

## Import BIN file
bin <- read_BIN2R(file = "D-CA1_100-150_2mm_a.BIN", fastForward = TRUE, 
                  txtProgressBar = FALSE, verbose = FALSE)

## Calculate equivalent doses
de <- analyse_SAR.CWOSL(object = bin, 
                        signal.integral.min = 1, signal.integral.max = 3, 
                        background.integral.min = 200, background.integral.max = 250, 
                        plot = FALSE, verbose = FALSE)

de_df <- get_RLum(de)
de_dist <- data.frame(de_df$De, de_df$De.Error)

## Minimum dose model
sigmab_range <- seq(from = 0.1, to = 0.2, by = 0.01)

mam <- lapply(sigmab_range, function(x) {
  calc_MinDose(de_dist, sigmab = x, log = TRUE, par = 3, plot = FALSE, verbose = FALSE)
})

mam_df <- do.call(rbind, get_RLum(mam))

## Results
mam_df
##          de   de_err ci_level ci_lower ci_upper par       sig        p0 mu
## 1  411.8003 38.39822     0.95 315.2927 465.8137   3 0.2736519 0.3571652 NA
## 2  416.3479 39.99034     0.95 316.3078 473.0700   3 0.2682456 0.3975128 NA
## 3  420.6191 41.60503     0.95 317.0109 480.1026   3 0.2617650 0.4325768 NA
## 4  424.8453 43.32688     0.95 317.4574 487.2988   3 0.2542803 0.4655925 NA
## 5  429.0137 45.31942     0.95 317.8340 495.4861   3 0.2455723 0.4963381 NA
## 6  433.1620 47.66035     0.95 318.1830 505.0116   3 0.2351025 0.5245625 NA
## 7  437.5172 49.80264     0.95 318.5311 513.7575   3 0.2222580 0.5524660 NA
## 8  442.4548 51.53531     0.95 318.9697 520.9881   3 0.2065044 0.5851851 NA
## 9  473.7277       NA     0.95       NA 531.1466   3 0.5430017 1.0000000 NA
## 10 473.8290       NA     0.95       NA 534.5057   3 0.5508769 1.0000000 NA
## 11 466.3760 54.96004     0.95 322.1196 537.5630   3 0.1347884 0.8500242 NA
##        Lmax      BIC
## 1  2.136488 7.709954
## 2  2.142909 7.697112
## 3  2.129751 7.723426
## 4  2.101564 7.779801
## 5  2.061718 7.859493
## 6  2.013184 7.956561
## 7  1.959375 8.064178
## 8  1.903985 8.174959
## 9  1.666726 8.649477
## 10 1.736842 8.509244
## 11 1.758732 8.465465
## Lowest BIC score
mam_df[which.min(mam_df$BIC), ]
##         de   de_err ci_level ci_lower ci_upper par       sig        p0 mu
## 2 416.3479 39.99034     0.95 316.3078   473.07   3 0.2682456 0.3975128 NA
##       Lmax      BIC
## 2 2.142909 7.697112