Package 'MBC' reference manual

Title:	Multivariate Bias Correction of Climate Model Outputs
Description:	Calibrate and apply multivariate bias correction algorithms for climate model simulations of multiple climate variables. Three methods described by Cannon (2016) <doi:10.1175/JCLI-D-15-0679.1> and Cannon (2018) <doi:10.1007/s00382-017-3580-6> are implemented — (i) MBC Pearson correlation (MBCp), (ii) MBC rank correlation (MBCr), and (iii) MBC N-dimensional PDF transform (MBCn) — as is the Rank Resampling for Distributions and Dependences (R2D2) method.
Authors:	Alex J. Cannon [aut, cre]
Maintainer:	Alex J. Cannon <[email protected]>
License:	GPL-2
Version:	0.10-7
Built:	2024-11-13 03:54:45 UTC
Source:	https://github.com/cran/MBC

Multivariate Bias Correction of Climate Model Outputs

Description

Calibrate and apply multivariate bias correction algorithms for climate model simulations of multiple climate variables. Three iterative methods are supported: (i) MBC Pearson correlation (MBCp), (ii) MBC rank correlation (MBCr), and (iii) MBC N-dimensional probability density function transform (MBCn).

The first two, MBCp and MBCr (Cannon, 2016), match marginal distributions and inter-variable dependence structure. Dependence structure can be measured either by the Pearson correlation (MBCp) or by the Spearman rank correlation (MBCr). The energy distance score (escore) is recommended for model selection.

The third, MBCn (Cannon, 2018), which operates on the full multivariate distribution, is more flexible and can be considered to be a multivariate analogue of univariate quantile mapping. All aspects of the observed distribution are transferred to the climate model simulations.

In each of the three methods, marginal distributions are corrected by the change-preserving quantile delta mapping (QDM) algorithm (Cannon et al., 2015). Finally, an implementation of the Rank Resampling for Distributions and Dependences (R2D2) method introduced by Vrac (2018) is also included.

An example application of the three MBC methods using the cccma dataset can be run via:

example(MBC, run.dontrun=TRUE)

Note: these functions apply bias correction to the supplied data without reference to other conditioning variables (e.g., time of year or season). The user must partition their data in a way that makes sense for a particular application. Furthermore, if MBCn is being applied to multiple spatial locations or grid cells, it is recommended that the same sequence of random rotations be used for each location. This can be achived by generating the random rotations first and then passing the sequence using the rot.seq argument:

rot.seq <- replicate(niterations, list(rot.random(nvars)))

bias_correction <- MBCn(..., rot.seq=rot.seq)

Details

Package:	MBC
Type:	Package
License:	GPL-2
LazyLoad:	yes

References

Cannon, A.J., 2018. Multivariate quantile mapping bias correction: An N-dimensional probability density function transform for climate model simulations of multiple variables. Climate Dynamics, 50(1-2):31-49. doi:10.1007/s00382-017-3580-6

Cannon, A.J., 2016. Multivariate bias correction of climate model output: Matching marginal distributions and inter-variable dependence structure. Journal of Climate, 29:7045-7064. doi:10.1175/JCLI-D-15-0679.1

Cannon, A.J., S.R. Sobie, and T.Q. Murdock, 2015. Bias correction of simulated precipitation by quantile mapping: How well do methods preserve relative changes in quantiles and extremes? Journal of Climate, 28:6938-6959. doi:10.1175/JCLI-D-14-00754.1

Francois, B., M. Vrac, A.J. Cannon, Y. Robin, and D. Allard, 2020. Multivariate bias corrections of climate simulations: Which benefits for which losses? Earth System Dynamics, 11:537-562. doi:10.5194/esd-11-537-2020

Vrac, M., 2018. Multivariate bias adjustment of high-dimensional climate simulations: the Rank Resampling for Distributions and Dependences (R2D2) bias correction. Hydrology and Earth System Sciences, 22:3175-3196. doi:10.5194/hess-22-3175-2018

Examples

## Not run: 
data(cccma)
set.seed(1)

# Univariate quantile mapping
qdm.c <- cccma$gcm.c*0
qdm.p <- cccma$gcm.p*0
for(i in seq(ncol(cccma$gcm.c))){
    fit.qdm <- QDM(o.c=cccma$rcm.c[,i], m.c=cccma$gcm.c[,i],
                   m.p=cccma$gcm.p[,i], ratio=cccma$ratio.seq[i],
                   trace=cccma$trace[i])
    qdm.c[,i] <- fit.qdm$mhat.c
    qdm.p[,i] <- fit.qdm$mhat.p
}

# Multivariate MBCp bias correction
fit.mbcp <- MBCp(o.c=cccma$rcm.c, m.c=cccma$gcm.c,
                 m.p=cccma$gcm.p, ratio.seq=cccma$ratio.seq,
                 trace=cccma$trace)
mbcp.c <- fit.mbcp$mhat.c
mbcp.p <- fit.mbcp$mhat.p

# Multivariate MBCr bias correction
fit.mbcr <- MBCr(o.c=cccma$rcm.c, m.c=cccma$gcm.c,
                 m.p=cccma$gcm.p, ratio.seq=cccma$ratio.seq,
                 trace=cccma$trace)
mbcr.c <- fit.mbcr$mhat.c
mbcr.p <- fit.mbcr$mhat.p

# Multivariate MBCn bias correction
fit.mbcn <- MBCn(o.c=cccma$rcm.c, m.c=cccma$gcm.c,
                 m.p=cccma$gcm.p, ratio.seq=cccma$ratio.seq,
                 trace=cccma$trace)
mbcn.c <- fit.mbcn$mhat.c
mbcn.p <- fit.mbcn$mhat.p
colnames(mbcn.c) <- colnames(mbcn.p) <-
    colnames(cccma$rcm.c)

# Correlation matrices (Pearson and Spearman)
# MBCp
dev.new()
par(mfrow=c(2, 2))
plot(c(cor(cccma$rcm.c)), c(cor(qdm.c)), col='black',
     pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCp',
     main='Pearson correlation\nMBCp calibration')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.c)), c(cor(mbcp.c)), col='red')
plot(c(cor(cccma$rcm.p)), c(cor(qdm.p)),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCp',
     main='Pearson correlation\nMBCp evaluation')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.p)), c(cor(mbcp.p)), col='red')
plot(c(cor(cccma$rcm.c, m='s')), c(cor(qdm.c, m='s')),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCp',
     main='Spearman correlation\nMBCp calibration')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.c, m='s')), c(cor(mbcp.c, m='s')),
       col='red')
plot(c(cor(cccma$rcm.p, m='s')), c(cor(qdm.p, m='s')),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCp',
     main='Spearman correlation\nMBCp evaluation')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.p, m='s')), c(cor(mbcp.p, m='s')),
       col='red')

# MBCr
dev.new()
par(mfrow=c(2, 2))
plot(c(cor(cccma$rcm.c)), c(cor(qdm.c)), col='black',
     pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCr',
     main='Pearson correlation\nMBCr calibration')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.c)), c(cor(mbcr.c)), col='blue')
plot(c(cor(cccma$rcm.p)), c(cor(qdm.p)),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCr',
     main='Pearson correlation\nMBCr evaluation')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.p)), c(cor(mbcr.p)), col='blue')
plot(c(cor(cccma$rcm.c, m='s')), c(cor(qdm.c, m='s')),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCr',
     main='Spearman correlation\nMBCr calibration')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.c, m='s')), c(cor(mbcr.c, m='s')),
       col='blue')
plot(c(cor(cccma$rcm.p, m='s')), c(cor(qdm.p, m='s')),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCr',
     main='Spearman correlation\nMBCr evaluation')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.p, m='s')), c(cor(mbcr.p, m='s')),
       col='blue')

# MBCn
dev.new()
par(mfrow=c(2, 2))
plot(c(cor(cccma$rcm.c)), c(cor(qdm.c)), col='black',
     pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCn',
     main='Pearson correlation\nMBCn calibration')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.c)), c(cor(mbcn.c)), col='orange')
plot(c(cor(cccma$rcm.p)), c(cor(qdm.p)),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCn',
     main='Pearson correlation\nMBCn evaluation')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.p)), c(cor(mbcn.p)), col='orange')
plot(c(cor(cccma$rcm.c, m='s')), c(cor(qdm.c, m='s')),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCn',
     main='Spearman correlation\nMBCn calibration')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.c, m='s')), c(cor(mbcn.c, m='s')),
       col='orange')
plot(c(cor(cccma$rcm.p, m='s')), c(cor(qdm.p, m='s')),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCn',
     main='Spearman correlation\nMBCn evaluation')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.p, m='s')), c(cor(mbcn.p, m='s')),
       col='orange')

# Pairwise scatterplots
dev.new()
pairs(cccma$gcm.c, main='CanESM2 calibration', col='#0000001A')
dev.new()
pairs(cccma$rcm.c, main='CanRCM4 calibration', col='#0000001A')
dev.new()
pairs(qdm.c, main='QDM calibration', col='#0000001A')
dev.new()
pairs(mbcp.c, main='MBCp calibration', col='#FF00001A')
dev.new()
pairs(mbcr.c, main='MBCr calibration', col='#0000FF1A')
dev.new()
pairs(mbcn.c, main='MBCn calibration', col='#FFA5001A')

# Energy distance skill score relative to univariate QDM
escore.qdm <- escore(cccma$rcm.p, qdm.p, scale.x=TRUE)
escore.mbcp <- escore(cccma$rcm.p, mbcp.p, scale.x=TRUE)
escore.mbcr <- escore(cccma$rcm.p, mbcr.p, scale.x=TRUE)
escore.mbcn <- escore(cccma$rcm.p, mbcn.p, scale.x=TRUE)

cat('ESS (MBCp):', 1-escore.mbcp/escore.qdm, '\n')
cat('ESS (MBCr):', 1-escore.mbcr/escore.qdm, '\n')
cat('ESS (MBCn):', 1-escore.mbcn/escore.qdm, '\n')

## End(Not run)
## Not run: 
data(cccma)
set.seed(1)

# Univariate quantile mapping
qdm.c <- cccma$gcm.c*0
qdm.p <- cccma$gcm.p*0
for(i in seq(ncol(cccma$gcm.c))){
    fit.qdm <- QDM(o.c=cccma$rcm.c[,i], m.c=cccma$gcm.c[,i],
                   m.p=cccma$gcm.p[,i], ratio=cccma$ratio.seq[i],
                   trace=cccma$trace[i])
    qdm.c[,i] <- fit.qdm$mhat.c
    qdm.p[,i] <- fit.qdm$mhat.p
}

# Multivariate MBCp bias correction
fit.mbcp <- MBCp(o.c=cccma$rcm.c, m.c=cccma$gcm.c,
                 m.p=cccma$gcm.p, ratio.seq=cccma$ratio.seq,
                 trace=cccma$trace)
mbcp.c <- fit.mbcp$mhat.c
mbcp.p <- fit.mbcp$mhat.p

# Multivariate MBCr bias correction
fit.mbcr <- MBCr(o.c=cccma$rcm.c, m.c=cccma$gcm.c,
                 m.p=cccma$gcm.p, ratio.seq=cccma$ratio.seq,
                 trace=cccma$trace)
mbcr.c <- fit.mbcr$mhat.c
mbcr.p <- fit.mbcr$mhat.p

# Multivariate MBCn bias correction
fit.mbcn <- MBCn(o.c=cccma$rcm.c, m.c=cccma$gcm.c,
                 m.p=cccma$gcm.p, ratio.seq=cccma$ratio.seq,
                 trace=cccma$trace)
mbcn.c <- fit.mbcn$mhat.c
mbcn.p <- fit.mbcn$mhat.p
colnames(mbcn.c) <- colnames(mbcn.p) <-
    colnames(cccma$rcm.c)

# Correlation matrices (Pearson and Spearman)
# MBCp
dev.new()
par(mfrow=c(2, 2))
plot(c(cor(cccma$rcm.c)), c(cor(qdm.c)), col='black',
     pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCp',
     main='Pearson correlation\nMBCp calibration')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.c)), c(cor(mbcp.c)), col='red')
plot(c(cor(cccma$rcm.p)), c(cor(qdm.p)),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCp',
     main='Pearson correlation\nMBCp evaluation')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.p)), c(cor(mbcp.p)), col='red')
plot(c(cor(cccma$rcm.c, m='s')), c(cor(qdm.c, m='s')),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCp',
     main='Spearman correlation\nMBCp calibration')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.c, m='s')), c(cor(mbcp.c, m='s')),
       col='red')
plot(c(cor(cccma$rcm.p, m='s')), c(cor(qdm.p, m='s')),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCp',
     main='Spearman correlation\nMBCp evaluation')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.p, m='s')), c(cor(mbcp.p, m='s')),
       col='red')

# MBCr
dev.new()
par(mfrow=c(2, 2))
plot(c(cor(cccma$rcm.c)), c(cor(qdm.c)), col='black',
     pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCr',
     main='Pearson correlation\nMBCr calibration')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.c)), c(cor(mbcr.c)), col='blue')
plot(c(cor(cccma$rcm.p)), c(cor(qdm.p)),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCr',
     main='Pearson correlation\nMBCr evaluation')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.p)), c(cor(mbcr.p)), col='blue')
plot(c(cor(cccma$rcm.c, m='s')), c(cor(qdm.c, m='s')),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCr',
     main='Spearman correlation\nMBCr calibration')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.c, m='s')), c(cor(mbcr.c, m='s')),
       col='blue')
plot(c(cor(cccma$rcm.p, m='s')), c(cor(qdm.p, m='s')),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCr',
     main='Spearman correlation\nMBCr evaluation')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.p, m='s')), c(cor(mbcr.p, m='s')),
       col='blue')

# MBCn
dev.new()
par(mfrow=c(2, 2))
plot(c(cor(cccma$rcm.c)), c(cor(qdm.c)), col='black',
     pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCn',
     main='Pearson correlation\nMBCn calibration')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.c)), c(cor(mbcn.c)), col='orange')
plot(c(cor(cccma$rcm.p)), c(cor(qdm.p)),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCn',
     main='Pearson correlation\nMBCn evaluation')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.p)), c(cor(mbcn.p)), col='orange')
plot(c(cor(cccma$rcm.c, m='s')), c(cor(qdm.c, m='s')),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCn',
     main='Spearman correlation\nMBCn calibration')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.c, m='s')), c(cor(mbcn.c, m='s')),
       col='orange')
plot(c(cor(cccma$rcm.p, m='s')), c(cor(qdm.p, m='s')),
     col='black', pch=19, xlim=c(-1, 1), ylim=c(-1, 1),
     xlab='CanRCM4', ylab='CanESM2 MBCn',
     main='Spearman correlation\nMBCn evaluation')
abline(0, 1)
grid()
points(c(cor(cccma$rcm.p, m='s')), c(cor(mbcn.p, m='s')),
       col='orange')

# Pairwise scatterplots
dev.new()
pairs(cccma$gcm.c, main='CanESM2 calibration', col='#0000001A')
dev.new()
pairs(cccma$rcm.c, main='CanRCM4 calibration', col='#0000001A')
dev.new()
pairs(qdm.c, main='QDM calibration', col='#0000001A')
dev.new()
pairs(mbcp.c, main='MBCp calibration', col='#FF00001A')
dev.new()
pairs(mbcr.c, main='MBCr calibration', col='#0000FF1A')
dev.new()
pairs(mbcn.c, main='MBCn calibration', col='#FFA5001A')

# Energy distance skill score relative to univariate QDM
escore.qdm <- escore(cccma$rcm.p, qdm.p, scale.x=TRUE)
escore.mbcp <- escore(cccma$rcm.p, mbcp.p, scale.x=TRUE)
escore.mbcr <- escore(cccma$rcm.p, mbcr.p, scale.x=TRUE)
escore.mbcn <- escore(cccma$rcm.p, mbcn.p, scale.x=TRUE)

cat('ESS (MBCp):', 1-escore.mbcp/escore.qdm, '\n')
cat('ESS (MBCr):', 1-escore.mbcr/escore.qdm, '\n')
cat('ESS (MBCn):', 1-escore.mbcn/escore.qdm, '\n')

## End(Not run)

Sample CanESM2 and CanRCM4 data

Description

Sample CanESM2 (T63 grid) and CanRCM4 (0.22-deg grid) data (122.5 deg W, 50 deg N).

pr: precipitation (mm day-1) 
tas: average surface temperature (deg. C)
dtr: diurnal temperature range (deg. C)
sfcWind: surface wind speed (m s-1)
ps: surface pressure (ps)
huss: surface specific humidity (kg kg-1)
rsds: surface downwelling shortwave radiation (W m-2)
rlds: surface downwelling longwave radiation (W m-2)

Value

a list of with elements consisting of:

`gcm.c`	matrix of CanESM2 variables for the calibration period.
`gcm.p`	matrix of CanESM2 variables for the validation period.
`rcm.c`	matrix of CanRCM4 variables for the calibration period.
`rcm.p`	matrix of CanRCM4 variables for the validation period.
`ratio.seq`	vector of logical values indicating if samples are of a ratio quantity.
`trace`	numeric values indicating trace thresholds for each ratio quantity.

Energy distance score

Description

Calculate the energy distance score measuring the statistical discrepancy between samples x and y from two multivariate distributions.

Usage

escore(x, y, scale.x=FALSE, n.cases=NULL, alpha=1, method='cluster')
escore(x, y, scale.x=FALSE, n.cases=NULL, alpha=1, method='cluster')

Arguments

`x`	numeric matrix.
`y`	numeric matrix.
`scale.x`	logical indicating whether data should be standardized based on `x`.
`n.cases`	the number of sub-sampled cases; `NULL` uses all data.
`alpha`	distance exponent in (0,2]
`method`	method used to weight the statistics

References

Székely, G.J. and M.L. Rizzo, 2004. Testing for equal distributions in high dimension, InterStat, November (5).

Székely, G.J. and M.L. Rizzo, 2013. Energy statistics: statistics based on distances. Journal of Statistical Planning and Inference, 143(8):1249-1272. doi:10.1016/j.jspi.2013.03.018

Rizzo, M.L. and G.L. Székely, 2016. Energy distance. Wiley Interdisciplinary Reviews: Computational Statistics, 8(1):27-38.

Multivariate bias correction (N-pdft)

Description

Multivariate bias correction that matches the multivariate distribution using QDM and the N-dimensional probability density function transform (N-pdft) following Cannon (2018).

Usage

MBCn(o.c, m.c, m.p, iter=30, ratio.seq=rep(FALSE, ncol(o.c)),
     trace=0.05, trace.calc=0.5*trace, jitter.factor=0,
     n.tau=NULL, ratio.max=2, ratio.max.trace=10*trace,
     ties='first', qmap.precalc=FALSE, rot.seq=NULL,
     silent=FALSE, n.escore=0, return.all=FALSE, subsample=NULL,
     pp.type=7)
MBCn(o.c, m.c, m.p, iter=30, ratio.seq=rep(FALSE, ncol(o.c)),
     trace=0.05, trace.calc=0.5*trace, jitter.factor=0,
     n.tau=NULL, ratio.max=2, ratio.max.trace=10*trace,
     ties='first', qmap.precalc=FALSE, rot.seq=NULL,
     silent=FALSE, n.escore=0, return.all=FALSE, subsample=NULL,
     pp.type=7)

Arguments

`o.c`	matrix of observed samples during the calibration period.
`m.c`	matrix of model outputs during the calibration period.
`m.p`	matrix of model outputs during the projected period.
`iter`	maximum number of algorithm iterations.
`ratio.seq`	vector of logical values indicating if samples are of a ratio quantity (e.g., precipitation).
`trace`	numeric values indicating thresholds below which values of a ratio quantity (e.g., `ratio=TRUE`) should be considered exact zeros.
`trace.calc`	numeric values of thresholds used internally when handling of exact zeros; defaults to one half of `trace`.
`jitter.factor`	optional strength of jittering to be applied when quantities are quantized.
`n.tau`	number of quantiles used in the quantile mapping; `NULL` equals the length of the `m.p` series.
`ratio.max`	numeric values indicating the maximum proportional changes allowed for ratio quantities below the `ratio.max.trace` threshold.
`ratio.max.trace`	numeric values of trace thresholds used to constrain the proportional change in ratio quantities to `ratio.max`; defaults to ten times `trace`.
`ties`	method used to handle ties when calculating ordinal ranks.
`qmap.precalc`	logical value indicating if `m.c` and `m.p` are outputs from `QDM`.
`rot.seq`	use a supplied list of random rotation matrices. `NULL` generates on the fly.
`silent`	logical value indicating if algorithm progress should be reported.
`n.escore`	number of cases used to calculate the energy distance when monitoring convergence.
`return.all`	logical value indicating whether results from all iterations are returned.
`subsample`	use `subsample` draws of size `n.tau` to calculate initial empirical quantiles; if `NULL`, calculate normally.
`pp.type`	type of plotting position used in `quantile`.

Value

a list of with elements consisting of:

`mhat.c`	matrix of bias corrected `m.c` values for the calibration period.
`mhat.p`	matrix of bias corrected `m.p` values for the projection period.

References

Pitié, F., A.C. Kokaram, and R. Dahyot, 2005. N-dimensional probability density function transfer and its application to color transfer. In Tenth IEEE International Conference on Computer Vision, 2005. ICCV 2005. (Vol. 2, pp. 1434-1439). IEEE.

Pitié, F., A.C. Kokaram, and R. Dahyot, 2007. Automated colour grading using colour distribution transfer. Computer Vision and Image Understanding, 107(1):123-137.

Multivariate bias correction (Pearson correlation)

Description

Multivariate bias correction that matches marginal distributions using QDM and the Pearson correlation dependence structure following Cannon (2016).

Usage

MBCp(o.c, m.c, m.p, iter=20, cor.thresh=1e-4,
     ratio.seq=rep(FALSE, ncol(o.c)), trace=0.05,
     trace.calc=0.5*trace, jitter.factor=0, n.tau=NULL,
     ratio.max=2, ratio.max.trace=10*trace, ties='first',
     qmap.precalc=FALSE, silent=FALSE, subsample=NULL,
     pp.type=7)
MBCp(o.c, m.c, m.p, iter=20, cor.thresh=1e-4,
     ratio.seq=rep(FALSE, ncol(o.c)), trace=0.05,
     trace.calc=0.5*trace, jitter.factor=0, n.tau=NULL,
     ratio.max=2, ratio.max.trace=10*trace, ties='first',
     qmap.precalc=FALSE, silent=FALSE, subsample=NULL,
     pp.type=7)

Arguments

`o.c`	matrix of observed samples during the calibration period.
`m.c`	matrix of model outputs during the calibration period.
`m.p`	matrix of model outputs during the projected period.
`iter`	maximum number of algorithm iterations.
`cor.thresh`	if greater than zero, a threshold indicating the change in magnitude of Pearson correlations required for convergence.
`ratio.seq`	vector of logical values indicating if samples are of a ratio quantity (e.g., precipitation).
`trace`	numeric values indicating thresholds below which values of a ratio quantity (e.g., `ratio=TRUE`) should be considered exact zeros.
`trace.calc`	numeric values of thresholds used internally when handling of exact zeros; defaults to one half of `trace`.
`jitter.factor`	optional strength of jittering to be applied when quantities are quantized.
`n.tau`	number of quantiles used in the quantile mapping; `NULL` equals the length of the `m.p` series.
`ratio.max`	numeric values indicating the maximum proportional changes allowed for ratio quantities below the `ratio.max.trace` threshold.
`ratio.max.trace`	numeric values of trace thresholds used to constrain the proportional change in ratio quantities to `ratio.max`; defaults to ten times `trace`.
`ties`	method used to handle ties when calculating ordinal ranks.
`qmap.precalc`	logical value indicating if `m.c` and `m.p` are outputs from `QDM`.
`silent`	logical value indicating if algorithm progress should be reported.
`subsample`	use `subsample` draws of size `n.tau` to calculate initial empirical quantiles; if `NULL`, calculate normally.
`pp.type`	type of plotting position used in `quantile`.

Value

a list of with elements consisting of:

`mhat.c`	matrix of bias corrected `m.c` values for the calibration period.
`mhat.p`	matrix of bias corrected `m.p` values for the projection period.

References

Multivariate bias correction (Spearman rank correlation)

Description

Multivariate bias correction that matches marginal distributions using QDM and the Spearman rank correlation dependence structure following Cannon (2016).

Usage

MBCr(o.c, m.c, m.p, iter=20, cor.thresh=1e-4,
     ratio.seq=rep(FALSE, ncol(o.c)), trace=0.05,
     trace.calc=0.5*trace, jitter.factor=0, n.tau=NULL,
     ratio.max=2, ratio.max.trace=10*trace, ties='first',
     qmap.precalc=FALSE, silent=FALSE, subsample=NULL,
     pp.type=7)
MBCr(o.c, m.c, m.p, iter=20, cor.thresh=1e-4,
     ratio.seq=rep(FALSE, ncol(o.c)), trace=0.05,
     trace.calc=0.5*trace, jitter.factor=0, n.tau=NULL,
     ratio.max=2, ratio.max.trace=10*trace, ties='first',
     qmap.precalc=FALSE, silent=FALSE, subsample=NULL,
     pp.type=7)

Arguments

`o.c`	matrix of observed samples during the calibration period.
`m.c`	matrix of model outputs during the calibration period.
`m.p`	matrix of model outputs during the projected period.
`iter`	maximum number of algorithm iterations.
`cor.thresh`	if greater than zero, a threshold indicating the change in magnitude of Spearman rank correlations required for convergence.
`ratio.seq`	vector of logical values indicating if samples are of a ratio quantity (e.g., precipitation).
`trace`	numeric values indicating thresholds below which values of a ratio quantity (e.g., `ratio=TRUE`) should be considered exact zeros.
`trace.calc`	numeric values of thresholds used internally when handling of exact zeros; defaults to one half of `trace`.
`jitter.factor`	optional strength of jittering to be applied when quantities are quantized.
`n.tau`	number of quantiles used in the quantile mapping; `NULL` equals the length of the `m.p` series.
`ratio.max`	numeric values indicating the maximum proportional changes allowed for ratio quantities below the `ratio.max.trace` threshold.
`ratio.max.trace`	numeric values of trace thresholds used to constrain the proportional change in ratio quantities to `ratio.max`; defaults to ten times `trace`.
`ties`	method used to handle ties when calculating ordinal ranks.
`qmap.precalc`	logical value indicating if `m.c` and `m.p` are outputs from `QDM`.
`silent`	logical value indicating if algorithm progress should be reported.
`subsample`	use `subsample` draws of size `n.tau` to calculate empirical quantiles; if `NULL`, calculate normally.
`pp.type`	type of plotting position used in `quantile`.

Value

a list of with elements consisting of:

`mhat.c`	matrix of bias corrected `m.c` values for the calibration period.
`mhat.p`	matrix of bias corrected `m.p` values for the projection period.

References

Multivariate linear rescaling using Cholesky decomposition

Description

Multivariate linear bias correction based on Cholesky decomposition of the covariance matrix following Scheuer and Stoller (1962) and Bürger et al. (2011). Bias correction matches the multivariate mean and covariance structure.

Usage

MRS(o.c, m.c, m.p, o.c.chol=NULL, o.p.chol=NULL, m.c.chol=NULL,
    m.p.chol=NULL)
MRS(o.c, m.c, m.p, o.c.chol=NULL, o.p.chol=NULL, m.c.chol=NULL,
    m.p.chol=NULL)

Arguments

`o.c`	matrix of observed samples during the calibration period.
`m.c`	matrix of model outputs during the calibration period.
`m.p`	matrix of model outputs during the projected period.
`o.c.chol`	precalculated Cholesky decomposition of the `o.c` covariance matrix; `NULL` calculates internally.
`o.p.chol`	precalculated Cholesky decomposition of the target `o.p` covariance matrix; `NULL` defaults to `o.c.chol`.
`m.c.chol`	precalculated Cholesky decomposition of the `m.c` covariance matrix; `NULL` calculates internally.
`m.p.chol`	precalculated Cholesky decomposition of the `m.p` covariance matrix; `NULL` calculates internally.

Value

a list of with elements consisting of:

`mhat.c`	matrix of bias corrected `m.c` values for the calibration period.
`mhat.p`	matrix of bias corrected `m.p` values for the projection period.

References

Scheuer, E.M. and D.S. Stoller, 1962. On the generation of normal random vectors. Technometrics, 4(2):278-281.

Bürger, G., J. Schulla, and A.T. Werner, 2011. Estimates of future flow, including extremes, of the Columbia River headwaters. Water Resources Research, 47(10):W10520. doi:10.1029/2010WR009716

Univariate bias correction via quantile delta mapping

Description

Univariate bias correction based on the quantile delta mapping QDM version of the quantile mapping algorithm from Cannon et al. (2015). QDM constrains model-projected changes in quantiles to be preserved following bias correction by quantile mapping.

Usage

QDM(o.c, m.c, m.p, ratio=FALSE, trace=0.05, trace.calc=0.5*trace,
    jitter.factor=0, n.tau=NULL, ratio.max=2,
    ratio.max.trace=10*trace, ECBC=FALSE, ties='first',
    subsample=NULL, pp.type=7)
QDM(o.c, m.c, m.p, ratio=FALSE, trace=0.05, trace.calc=0.5*trace,
    jitter.factor=0, n.tau=NULL, ratio.max=2,
    ratio.max.trace=10*trace, ECBC=FALSE, ties='first',
    subsample=NULL, pp.type=7)

Arguments

`o.c`	vector of observed samples during the calibration period.
`m.c`	vector of model outputs during the calibration period.
`m.p`	vector of model outputs during the projected period.
`ratio`	logical value indicating if samples are of a ratio quantity (e.g., precipitation).
`trace`	numeric value indicating the threshold below which values of a ratio quantity (e.g., `ratio=TRUE`) should be considered exact zeros.
`trace.calc`	numeric value of a threshold used internally when handling of exact zeros; defaults to one half of `trace`.
`jitter.factor`	optional strength of jittering to be applied when quantities are quantized.
`n.tau`	number of quantiles used in the quantile mapping; `NULL` equals the length of the `m.p` series.
`ratio.max`	numeric value indicating the maximum proportional change allowed for ratio quantities below the `ratio.max.trace` threshold.
`ratio.max.trace`	numeric value of a trace threshold used to constrain the proportional change in ratio quantities to `ratio.max`; defaults to ten times `trace`.
`ECBC`	logical value indicating whether `mhat.p` outputs should be ordered according to `o.c` ranks, i.e., as in the empirical copula-bias correction (ECBC) algorithm.
`ties`	method used to handle ties when calculating ordinal ranks.
`subsample`	use `subsample` draws of size `n.tau` to calculate empirical quantiles; if `NULL`, calculate normally.
`pp.type`	type of plotting position used in `quantile`.

Value

a list of with elements consisting of:

`mhat.c`	vector of bias corrected `m.c` values for the calibration period.
`mhat.p`	vector of bias corrected `m.p` values for the projection period.

References

Multivariate bias correction (R2D2)

Description

Multivariate bias correction that matches the multivariate distribution using QDM and the R2D2 algorithm following Vrac (2018).

Usage

R2D2(o.c, m.c, m.p, ref.column=1, ratio.seq=rep(FALSE, ncol(o.c)),
     trace=0.05, trace.calc=0.5*trace, jitter.factor=0,
     n.tau=NULL, ratio.max=2, ratio.max.trace=10*trace,
     ties='first', qmap.precalc=FALSE, subsample=NULL, pp.type=7)
R2D2(o.c, m.c, m.p, ref.column=1, ratio.seq=rep(FALSE, ncol(o.c)),
     trace=0.05, trace.calc=0.5*trace, jitter.factor=0,
     n.tau=NULL, ratio.max=2, ratio.max.trace=10*trace,
     ties='first', qmap.precalc=FALSE, subsample=NULL, pp.type=7)

Arguments

`o.c`	matrix of observed samples during the calibration period.
`m.c`	matrix of model outputs during the calibration period.
`m.p`	matrix of model outputs during the projected period.
`ref.column`	index of the reference column used for the 1D nearest neighbour matching
`ratio.seq`	vector of logical values indicating if samples are of a ratio quantity (e.g., precipitation).
`trace`	numeric values indicating thresholds below which values of a ratio quantity (e.g., `ratio=TRUE`) should be considered exact zeros.
`trace.calc`	numeric values of thresholds used internally when handling of exact zeros; defaults to one half of `trace`.
`jitter.factor`	optional strength of jittering to be applied when quantities are quantized.
`n.tau`	number of quantiles used in the quantile mapping; `NULL` equals the length of the `m.p` series.
`ratio.max`	numeric values indicating the maximum proportional changes allowed for ratio quantities below the `ratio.max.trace` threshold.
`ratio.max.trace`	numeric values of trace thresholds used to constrain the proportional change in ratio quantities to `ratio.max`; defaults to ten times `trace`.
`ties`	method used to handle ties when calculating ordinal ranks.
`qmap.precalc`	logical value indicating if `m.c` and `m.p` are outputs from `QDM`.
`subsample`	use `subsample` draws of size `n.tau` to calculate initial empirical quantiles; if `NULL`, calculate normally.
`pp.type`	type of plotting position used in `quantile`.

Value

a list of with elements consisting of:

`mhat.c`	matrix of bias corrected `m.c` values for the calibration period.
`mhat.p`	matrix of bias corrected `m.p` values for the projection period.

References

Random orthogonal rotation

Description

Generate a k-dimensional random orthogonal rotation matrix.

Usage

rot.random(k)
rot.random(k)

Arguments

`k`	the number of dimensions.

References

Mezzadri, F. 2007. How to generate random matrices from the classical compact groups, Notices of the American Mathematical Society, 54:592–604.

Package 'MBC'

Help Index

Multivariate Bias Correction of Climate Model Outputs

Description

Details

References

See Also

Examples

Sample CanESM2 and CanRCM4 data

Description

Value

Energy distance score

Description

Usage

Arguments

References

Multivariate bias correction (N-pdft)

Description

Usage

Arguments

Value

References

See Also

Multivariate bias correction (Pearson correlation)

Description

Usage

Arguments

Value

References

See Also

Multivariate bias correction (Spearman rank correlation)

Description

Usage

Arguments

Value

References

See Also

Multivariate linear rescaling using Cholesky decomposition

Description

Usage

Arguments

Value

References

See Also

Univariate bias correction via quantile delta mapping

Description

Usage

Arguments

Value

References

See Also

Multivariate bias correction (R2D2)

Description

Usage

Arguments

Value

References

See Also

Random orthogonal rotation

Description

Usage

Arguments

References