Apply a SMITE model (m) to a 'forward' coral variable matrix (A) to produce predicted values for a reconstruction target (bhat).

SMITE.recon(
  A,
  m,
  Amu,
  Asd,
  bmu,
  bsd
)

Arguments

A

The 'forward' coral variable matrix (t x p; time in rows, coral parameters in columns).

m

A vector of SMITE model parameters (length p), ideally created from the SMITE.calib() function.

Amu

A vector of length p giving the mean values of each coral variable in the calibration dataset.

Asd

A vector of length p giving the standard deviations of each coral variable in the calibration dataset.

bmu

A single-element atomic vector giving the mean of the reconstruction target in the calibration dataset.

bsd

A single-element atomic vector giving the standard deviation of the reconstruction target in the calibration dataset.

Value

Returns a vector of predicted reconstruction target values (bhat). Minimum errors (i.e., not included age model errors) for bhat can be estimated from the SEP calculated in the SMITE.calib() function.

References

Hughes, H. P., Thompson, D. M., Foster, G. L., Lees, J., Surge, D. (in revision). “Synthetic and Practical Reconstructions of SST and pH Using the Novel Multiproxy SMITE Method.” PLOS ONE.

Author

Hunter P. Hughes

Examples

# Load data from Hughes et al. (in revision)
library(SMITER)
data(BMDA_1B_Comp)
data(BMDA_1B_EComp)

# Designate forward matrix (A), reconstruction target (b), and corresponding errors (Ae, be).

proxy_names <- names(bmda_1b_comp[,c(5:ncol(bmda_1b_comp))])

A <- bmda_1b_comp[,proxy_names]
Ae <- bmda_1b_ecomp[,proxy_names]
b <- bmda_1b_comp[,'Temp']
be <- rep(0.02, length(b))

# Execute SMITE Calibration with no truncation (i.e., all information retained).

result <- SMITE.calib(A = A, b = b, Ae = Ae, be = be, eigenclean = ncol(A))

# Extract SMITE model parameters #
m <- result$m$m.mu

# Reconstruction (should be nearly identical to result$recon$bhat.mu, only minor differences due to bootstrapping)
bhat <- SMITE.recon(
  A = A,
  m = m,
  Amu = colMeans(A),
  Asd = apply(A, 2, function(x) sd(x)),
  bmu = mean(b),
  bsd = sd(b)
)

cbind(bhat, result$recon$bhat.mu)
#>           [,1]     [,2]
#>  [1,] 29.08020 28.86139
#>  [2,] 26.21584 26.49063
#>  [3,] 26.71475 26.53376
#>  [4,] 25.59253 25.40753
#>  [5,] 23.70693 23.54265
#>  [6,] 21.53489 21.46609
#>  [7,] 20.86025 20.94818
#>  [8,] 20.28968 20.51567
#>  [9,] 22.62603 22.59325
#> [10,] 25.28751 24.86350
#> [11,] 27.11015 26.85822
#> [12,] 28.28851 27.97278
#> [13,] 28.79140 28.35513
#> [14,] 27.62780 27.59611
#> [15,] 26.79474 26.38086
#> [16,] 24.75985 24.90373
#> [17,] 21.91731 22.23570
#> [18,] 19.76879 20.33556
#> [19,] 19.89459 20.44850
#> [20,] 20.92776 21.22034
#> [21,] 22.71602 22.57487
#> [22,] 21.74459 21.89511
#> [23,] 23.98782 23.90824
#> [24,] 27.82003 27.44491
#> [25,] 28.49472 28.26452
#> [26,] 28.10520 27.74512
#> [27,] 25.92628 25.78140
#> [28,] 23.20079 23.35733
#> [29,] 21.68172 21.96019
#> [30,] 19.71963 20.16831
#> [31,] 19.43738 19.79337
#> [32,] 20.09220 20.16216
#> [33,] 20.95260 21.08338