ctmva (version 1.6.0)

ccor_posim: Credible interval estimation for curve correlation based on posterior simulation

Description

Inputs raw data representing two curves, and computes a credible interval for the curve correlation between them simulating from the approximate posterior distribution of the joint spline coefficient vector.

Usage

ccor_posim(
  y,
  time,
  curve = NULL,
  method,
  k = 15,
  conf = 0.95,
  ndraw = 999,
  min.overlap = 0
)

Value

A list with components

cor

curve correlation

model

the model for the two curves (if method=="mvn"), or a list of the two curves' models (if method=="indep")

bsb

B-spline basis (from package fda) used for the curves

Vc.fda

corrected posterior covariance matrix for the coefficients with respect to the B-spline basis bsb (see the component $Vc in gamObject)

sims

curve correlations for the ndraw draws from the posterior distribution

ci

credible interval for the curve correlation

Arguments

y, time, curve, k, min.overlap

see ccor

method

"indep" (curves fitted independently) or "mvn" (curves fitted together, with error correlation estimated based on multivariate normal likelihood)

conf

confidence level

ndraw

number of draws from posterior distribution of spline coefficient vector

Author

Philip Tzvi Reiss <[e-mail]>, Noemi Foa, Dror Arbiv and Biplab Paul <[e-mail]>

See Also

ccor, ccor_boot, mvn

Examples

## Not run:
if (interactive () &&
    requireNamespace("wbwdi", quietly = TRUE) &&
    requireNamespace("ggplot2", quietly = TRUE) &&
    requireNamespace("dplyr", quietly = TRUE)) {

    # Curve correlation of per capita GDP and fertility rate in Paraguay
    wdi_dat <- wbwdi::wdi_get(entities = c("PRY"), start_year=1960, end_year=2023,
                       indicators = c("NY.GDP.PCAP.KD","SP.DYN.TFRT.IN"), format="wide") |>
        dplyr::rename(percapitaGDP = NY.GDP.PCAP.KD, fertility = SP.DYN.TFRT.IN)
    ggplot2::ggplot(wdi_dat, aes(percapitaGDP, fertility, color=year)) + geom_point()

    y <- as.matrix(wdi_dat[ , c("percapitaGDP", "fertility")])

    set.seed(345)

    ci <- list()
    ci[[1]] <- ccor_posim(y=y, time=wdi_dat$year, method="indep")
    ci[[2]] <- ccor_posim(y=y, time=wdi_dat$year, method="mvn")
    ci[[3]] <- ccor_boot(y=y, time=wdi_dat$year, ndraw=399)

    tabl <- matrix(NA, 3, 3)
    for (k in 1:3) tabl[k, ] <- c(ci[[k]]$cor, ci[[k]]$ci)
    dimnames(tabl) <- list(c("Posim_indep", "Posim_MVN", "Bootstrap"), c("Est","Lower95","Upper95"))
    round(tabl, 4)
}
## End(Not run)