From 3036e9555fcfd6228dc60153db95b4560dcd58f3 Mon Sep 17 00:00:00 2001 From: "J. Derek Tucker" Date: Fri, 17 Nov 2023 20:47:35 -0700 Subject: [PATCH] bugfixes --- R/elastic_changepoint.R | 593 ++++++++++++++++++----------------- man/elastic_amp_change_ff.Rd | 2 +- man/elastic_change_fpca.Rd | 2 +- man/elastic_ph_change_ff.Rd | 2 +- 4 files changed, 300 insertions(+), 299 deletions(-) diff --git a/R/elastic_changepoint.R b/R/elastic_changepoint.R index c5e71bb..85307fc 100644 --- a/R/elastic_changepoint.R +++ b/R/elastic_changepoint.R @@ -22,7 +22,7 @@ #' \item{WarpingMeanBefore}{mean warping function before changepoint} #' \item{WarpingMeanAfter}{mean warping function after changepoint} #' \item{change_fun}{amplitude change function} -#' \item{Sn}{test statisitc values} +#' \item{Sn}{test statistic values} #' \item{mu}{mean srsfs} #' \item{mu_f}{mean functions} #' @keywords srvf alignment changepoint @@ -32,107 +32,107 @@ #' @export elastic_amp_change_ff <- function(f, time, d = 1000, h = 0, smooth_data=FALSE, sparam=25, showplot = TRUE) { - if (smooth_data) f <- smooth.data(f, sparam) - - out <- time_warping(f, time, parallel = T, showplot=F) - - N1 <- dim(f)[2] - M <- dim(f)[1] - mu <- matrix(0, nrow = M, ncol = N1) - mu_f <- matrix(0, nrow = M, ncol = N1) - - # Compute Karcher mean for first i+1 functions - mu[, 1] <- out$qn[,1] - mu_f[, 1] <- out$fn[,1] - for (i in 2:N1) { - mu[, i] <- rowSums(out$qn[,1:i]) - mu_f[, i] <- rowSums(out$fn[,1:i]) - } - - # compute test statistic - Sn <- (1:N1) - Sn[1] <- 0 - for (j in (2:N1)) { - Sn[j] <- 1/M * pvecnorm(1/sqrt(N1)*(mu_f[,j-1]-(j/N1)*mu_f[,N1]))^2 + if (smooth_data) f <- smooth.data(f, sparam) + + out <- time_warping(f, time, parallel = T) + + N1 <- dim(f)[2] + M <- dim(f)[1] + mu <- matrix(0, nrow = M, ncol = N1) + mu_f <- matrix(0, nrow = M, ncol = N1) + + # Compute Karcher mean for first i+1 functions + mu[, 1] <- out$qn[,1] + mu_f[, 1] <- out$fn[,1] + for (i in 2:N1) { + mu[, i] <- rowSums(out$qn[,1:i]) + mu_f[, i] <- rowSums(out$fn[,1:i]) + } + + # compute test statistic + Sn <- (1:N1) + Sn[1] <- 0 + for (j in (2:N1)) { + Sn[j] <- 1/M * pvecnorm(1/sqrt(N1)*(mu_f[,j-1]-(j/N1)*mu_f[,N1]))^2 + } + k.star <- min(which(Sn == max(Sn))) + Tn <- max(Sn) + + # compute means on either side of the changepoint + dat.a <- f[, 1:k.star] + dat.b <- f[, (k.star + 1):N1] + warp.a <- out$warping_functions[,1:k.star] + warp.b <- out$warping_functions[,(k.star + 1):N1] + mean.a <- rowMeans(out$fn[,1:k.star]) + mean.b <- rowMeans(out$fn[,(k.star + 1):N1]) + out.gam.mean <- SqrtMean(warp.a) + warp_mean_a = out.gam.mean$gam_mu + out.gam.mean <- SqrtMean(warp.b) + warp_mean_b = out.gam.mean$gam_mu + delta <- mean.a - mean.b + + # center your data + centered_data = matrix(0, M, N1) + for (i in (1:(k.star-1))){ + centered_data[,i] = out$fn[,i] - mean.a + } + for (i in (k.star:N1)){ + centered_data[,i] = out$fn[,i] - mean.b + } + + D_mat <- LongRunCovMatrix(centered_data, h = h) + eigen_struct <- eigen(D_mat, symmetric = TRUE) + lambda <- Re(eigen_struct$values)/M + + asymp <- function(N) { + BridgeLam <- matrix(0, M, N) + for (j in (1:M)) { + BridgeLam[j, ] <- lambda[j] * (BBridge(0, 0, 0, 1, N - 1)^2) } - k.star <- min(which(Sn == max(Sn))) - Tn <- max(Sn) - - # compute means on either side of the changepoint - dat.a <- f[, 1:k.star] - dat.b <- f[, (k.star + 1):N1] - warp.a <- out$warping_functions[,1:k.star] - warp.b <- out$warping_functions[,(k.star + 1):N1] - mean.a <- rowMeans(out$fn[,1:k.star]) - mean.b <- rowMeans(out$fn[,(k.star + 1):N1]) - out.gam.mean <- SqrtMean(warp.a) - warp_mean_a = out.gam.mean$gam_mu - out.gam.mean <- SqrtMean(warp.b) - warp_mean_b = out.gam.mean$gam_mu - delta <- mean.a - mean.b - - # center your data - centered_data = matrix(0, M, N1) - for (i in (1:(k.star-1))){ - centerd_data[,i] = out$fn[,i] - mean.a + max(colSums(BridgeLam)) + } + + Values <- sapply(1:d, function(k) asymp(N1)) + z <- Tn <= Values + p <- length(z[z == TRUE]) / length(z) + + # Plot + if (showplot == TRUE) { + par(mfrow = c(1, 3)) + matplot(f, type = "l", col = "grey", main = "Functional Data", ylab = "values") + for (i in 1:ncol(dat.a)) { + lines(dat.a[, i], col = "pink") } - for (i in (k.star:N1)){ - centerd_data[,i] = out$fn[,i] - mean.b + for (i in 1:ncol(dat.b)) { + lines(dat.b[, i], col = "lightblue") } - - D_mat <- LongRunCovMatrix(centered_data, h = h) - eigen_struct <- eigen(D_mat, symmetric = TRUE) - lambda <- Re(eigen_struct$values)/M - - asymp <- function(N) { - BridgeLam <- matrix(0, M, N) - for (j in (1:M)) { - BridgeLam[j, ] <- lambda[j] * (BBridge(0, 0, 0, 1, N - 1)^2) - } - max(colSums(BridgeLam)) + lines(mean.b, col = "blue") + lines(mean.a, col = "red") + legend("topleft", c("before", "after"), col = c("blue", "red"), lty = c(1, 1), cex = 0.5) + plot(delta, type = "l", main = "Estimated Change Function", ylab = "values") + + matplot(out$warping_functions, type = "l", col = "grey", main = "Warping Functions", ylab = "values") + for (i in 1:ncol(warp.b)) { + lines(warp.b[, i], col = "pink") } - - Values <- sapply(1:d, function(k) asymp(N1)) - z <- Tn <= Values - p <- length(z[z == TRUE]) / length(z) - - # Plot - if (showplot == TRUE) { - par(mfrow = c(1, 3)) - matplot(f, type = "l", col = "grey", main = "Functional Data", ylab = "values") - for (i in 1:ncol(dat.a)) { - lines(dat.a[, i], col = "pink") - } - for (i in 1:ncol(dat.b)) { - lines(dat.b[, i], col = "lightblue") - } - lines(mean.b, col = "blue") - lines(mean.a, col = "red") - legend("topleft", c("before", "after"), col = c("blue", "red"), lty = c(1, 1), cex = 0.5) - plot(delta, type = "l", main = "Estimated Change Function", ylab = "values") - - matplot(out$warping_functions, type = "l", col = "grey", main = "Warping Functions", ylab = "values") - for (i in 1:ncol(warp.b)) { - lines(warp.b[, i], col = "pink") - } - for (i in 1:ncol(warp.a)) { - lines(warp.a[, i], col = "lightblue") - } - lines(warp_mean_a, col = "blue") - lines(warp_mean_b, col = "red") - legend("topleft", c("before", "after"), col = c("blue", "red"), lty = c(1, 1), cex = 0.5) + for (i in 1:ncol(warp.a)) { + lines(warp.a[, i], col = "lightblue") } - - out <- list( - pvalue = p, change = k.star, - DataBefore = dat.a, DataAfter = dat.b, - MeanBefore = mean.a, MeanAfter = mean.b, - WarpingMeanBefore = warp_mean_a, WarpingMeanAfter = warp_mean_b, - WarpingBefore = warp.a, WarpingAfter = warp.b, - change_fun = delta, Sn = Sn, mu=mu, mu_f=mu_f - ) - - return(out) + lines(warp_mean_a, col = "blue") + lines(warp_mean_b, col = "red") + legend("topleft", c("before", "after"), col = c("blue", "red"), lty = c(1, 1), cex = 0.5) + } + + out <- list( + pvalue = p, change = k.star, + DataBefore = dat.a, DataAfter = dat.b, + MeanBefore = mean.a, MeanAfter = mean.b, + WarpingMeanBefore = warp_mean_a, WarpingMeanAfter = warp_mean_b, + WarpingBefore = warp.a, WarpingAfter = warp.b, + change_fun = delta, Sn = Sn, mu=mu, mu_f=mu_f + ) + + return(out) } @@ -160,7 +160,7 @@ elastic_amp_change_ff <- function(f, time, d = 1000, h = 0, smooth_data=FALSE, s #' \item{WarpingMeanBefore}{mean warping function before changepoint} #' \item{WarpingMeanAfter}{mean warping function after changepoint} #' \item{change_fun}{amplitude change function} -#' \item{Sn}{test statisitc values} +#' \item{Sn}{test statistic values} #' \item{mu}{mean shooting vectors} #' @keywords srvf alignment changepoint #' @references J. D. Tucker and D. Yarger, “Elastic Functional Changepoint @@ -169,111 +169,111 @@ elastic_amp_change_ff <- function(f, time, d = 1000, h = 0, smooth_data=FALSE, s #' @export elastic_ph_change_ff <- function(f, time, d = 1000, h = 0, smooth_data=FALSE, sparam=25, showplot = TRUE) { - if (smooth_data) f <- smooth.data(f, sparam) - - out <- time_warping(f, time, parallel = T, showplot=F) - - N1 <- dim(f)[2] - M <- dim(f)[1] - mu <- matrix(0, nrow = M, ncol = N1) - - # Compute Karcher mean of warping functions - out = SqrtMean(out.warping_functions) - vec = out$vec - - mu[, 1] <- vec[,1] - for (i in 2:N) { - mu[, i] <- rowSums(vec[,1:i]) + if (smooth_data) f <- smooth.data(f, sparam) + + out <- time_warping(f, time, parallel = T) + + N1 <- dim(f)[2] + M <- dim(f)[1] + mu <- matrix(0, nrow = M, ncol = N1) + + # Compute Karcher mean of warping functions + out.mean = SqrtMean(out$warping_functions) + vec = out.mean$vec + + mu[, 1] <- vec[,1] + for (i in 2:N1) { + mu[, i] <- rowSums(vec[,1:i]) + } + + # compute test statistic + Sn <- (1:N1) + Sn[1] <- 0 + for (j in (2:N1)) { + Sn[j] <- 1/M * pvecnorm(1/sqrt(N1)*(mu[,j-1]-(j/N1)*mu[,N1]))^2 + } + k.star <- min(which(Sn == max(Sn))) + Tn <- max(Sn) + + # compute means on either side of the changepoint + dat.a <- f[, 1:k.star] + dat.b <- f[, (k.star + 1):N1] + warp.a <- out$warping_functions[,1:k.star] + warp.b <- out$warping_functions[,(k.star + 1):N1] + mean.a <- rowMeans(out$fn[,1:k.star]) + mean.b <- rowMeans(out$fn[,(k.star + 1):N1]) + out.gam.mean <- SqrtMean(warp.a) + warp_mean_a = out.gam.mean$gam_mu + mu.a = rowMeans(out.gam.mean$vec) + out.gam.mean <- SqrtMean(warp.b) + warp_mean_b = out.gam.mean$gam_mu + mu.b = rowMeans(out.gam.mean$vec) + delta <- mean.a - mean.b + + # center your data + centered_data = matrix(0, M, N1) + for (i in (1:(k.star-1))){ + centered_data[,i] = vec[,i] - mu.a + } + for (i in (k.star:N1)){ + centered_data[,i] = vec[,i] - mu.b + } + + # estimate eigenvalues of covariance operator + D_mat <- LongRunCovMatrix(centered_data, h = h) + eigen_struct <- eigen(D_mat, symmetric = TRUE) + lambda <- Re(eigen_struct$values)/M + + asymp <- function(N) { + BridgeLam <- matrix(0, M, N) + for (j in (1:M)) { + BridgeLam[j, ] <- lambda[j] * (BBridge(0, 0, 0, 1, N - 1)^2) } - - # compute test statistic - Sn <- (1:N1) - Sn[1] <- 0 - for (j in (2:N1)) { - Sn[j] <- 1/M * pvecnorm(1/sqrt(N1)*(mu[,j-1]-(j/N1)*mu[,N1]))^2 + max(colSums(BridgeLam)) + } + + Values <- sapply(1:M, function(k) asymp(N1)) + z <- Tn <= Values + p <- length(z[z == TRUE]) / length(z) + + # Plot + delta <- mean.a - mean.b + if (showplot == TRUE) { + par(mfrow = c(1, 3)) + matplot(f, type = "l", col = "grey", main = "Functional Data", ylab = "values") + for (i in 1:ncol(dat.a)) { + lines(dat.a[, i], col = "pink") } - k.star <- min(which(Sn == max(Sn))) - Tn <- max(Sn) - - # compute means on either side of the changepoint - dat.a <- f[, 1:k.star] - dat.b <- f[, (k.star + 1):N1] - warp.a <- out$warping_functions[,1:k.star] - warp.b <- out$warping_functions[,(k.star + 1):N1] - mean.a <- rowMeans(out$fn[,1:k.star]) - mean.b <- rowMeans(out$fn[,(k.star + 1):N1]) - out.gam.mean <- SqrtMean(warp.a) - warp_mean_a = out.gam.mean$gam_mu - mu.a = rowMeans(out.gam.mean$vec) - out.gam.mean <- SqrtMean(warp.b) - warp_mean_b = out.gam.mean$gam_mu - mu.b = rowMeans(out.gam.mean$vec) - delta <- mean.a - mean.b - - # center your data - centered_data = matrix(0, M, N1) - for (i in (1:(k.star-1))){ - centerd_data[,i] = vec[,i] - mu.a + for (i in 1:ncol(dat.b)) { + lines(dat.b[, i], col = "lightblue") } - for (i in (k.star:N1)){ - centerd_data[,i] = vec[,i] - mu.b + lines(mean.b, col = "blue") + lines(mean.a, col = "red") + legend("topleft", c("before", "after"), col = c("blue", "red"), lty = c(1, 1), cex = 0.5) + plot(delta, type = "l", main = "Estimated Change Function", ylab = "values") + + matplot(out$warping_functions, type = "l", col = "grey", main = "Warping Functions", ylab = "values") + for (i in 1:ncol(warp.b)) { + lines(warp.b[, i], col = "pink") } - - # estimate eigenvalues of covariance operator - D_mat <- LongRunCovMatrix(centered_data, h = h) - eigen_struct <- eigen(D_mat, symmetric = TRUE) - lambda <- Re(eigen_struct$values)/M - - asymp <- function(N) { - BridgeLam <- matrix(0, M, N) - for (j in (1:M)) { - BridgeLam[j, ] <- lambda[j] * (BBridge(0, 0, 0, 1, N - 1)^2) - } - max(colSums(BridgeLam)) + for (i in 1:ncol(warp.a)) { + lines(warp.a[, i], col = "lightblue") } - - Values <- sapply(1:M1, function(k) asymp(N)) - z <- Tn <= Values - p <- length(z[z == TRUE]) / length(z) - - # Plot - delta <- mean.a - mean.b - if (showplot == TRUE) { - par(mfrow = c(1, 3)) - matplot(f, type = "l", col = "grey", main = "Functional Data", ylab = "values") - for (i in 1:ncol(dat.a)) { - lines(dat.a[, i], col = "pink") - } - for (i in 1:ncol(dat.b)) { - lines(dat.b[, i], col = "lightblue") - } - lines(mean.b, col = "blue") - lines(mean.a, col = "red") - legend("topleft", c("before", "after"), col = c("blue", "red"), lty = c(1, 1), cex = 0.5) - plot(delta, type = "l", main = "Estimated Change Function", ylab = "values") - - matplot(out$warping_functions, type = "l", col = "grey", main = "Warping Functions", ylab = "values") - for (i in 1:ncol(warp.b)) { - lines(warp.b[, i], col = "pink") - } - for (i in 1:ncol(warp.a)) { - lines(warp.a[, i], col = "lightblue") - } - lines(warp_mean_a, col = "blue") - lines(warp_mean_b, col = "red") - legend("topleft", c("before", "after"), col = c("blue", "red"), lty = c(1, 1), cex = 0.5) - } - - out <- list( - pvalue = p, change = k.star, - DataBefore = dat.a, DataAfter = dat.b, - MeanBefore = mean.a, MeanAfter = mean.b, - WarpingMeanBefore = warp_mean_a, WarpingMeanAfter = warp_mean_b, - WarpingBefore = warp.a, WarpingAfter = warp.b, - change_fun = delta, Sn = Sn, mu=mu - ) - - return(out) + lines(warp_mean_a, col = "blue") + lines(warp_mean_b, col = "red") + legend("topleft", c("before", "after"), col = c("blue", "red"), lty = c(1, 1), cex = 0.5) + } + + out <- list( + pvalue = p, change = k.star, + DataBefore = dat.a, DataAfter = dat.b, + MeanBefore = mean.a, MeanAfter = mean.b, + WarpingMeanBefore = warp_mean_a, WarpingMeanAfter = warp_mean_b, + WarpingBefore = warp.a, WarpingAfter = warp.b, + change_fun = delta, Sn = Sn, mu=mu + ) + + return(out) } @@ -303,113 +303,114 @@ elastic_ph_change_ff <- function(f, time, d = 1000, h = 0, smooth_data=FALSE, sp #' \item{WarpingMeanBefore}{mean warping function before changepoint} #' \item{WarpingMeanAfter}{mean warping function after changepoint} #' \item{change_fun}{amplitude change function} -#' \item{Sn}{test statisitc values} +#' \item{Sn}{test statistic values} #' @keywords srvf alignment changepoint #' @references J. D. Tucker and D. Yarger, “Elastic Functional Changepoint #' Detection of Climate Impacts from Localized Sources”, Envirometrics, #' 10.1002/env.2826, 2023. #' @export elastic_change_fpca <- function(f, time, pca.method = "combined", pc = 0.95, d = 1000, n_pcs = 5, smooth_data=FALSE, sparam=25, showplot = TRUE) { - N1 = length(time) - if (smooth_data) f <- smooth.data(f, sparam) - - pca.method1 <- pca.method - pca.method <- pmatch(pca.method, c("combined", "vert", "horiz")) # 1 - combined, 2 - vert, 3 - horiz - if (is.na(pca.method)) { - stop("invalid method selection") + M = length(time) + N1 = dim(f)[2] + if (smooth_data) f <- smooth.data(f, sparam) + + pca.method1 <- pca.method + pca.method <- pmatch(pca.method, c("combined", "vert", "horiz")) # 1 - combined, 2 - vert, 3 - horiz + if (is.na(pca.method)) { + stop("invalid method selection") + } + + out <- time_warping(f, time, parallel = T) + + # Calculate PCA ----------------------------------------------------------- + if (pca.method == 1) + out.pca <- jointFPCA(out, M, showplot=F, C=NULL) + if (pca.method == 2) + out.pca <- vertFPCA(out, M, showplot=F) + if (pca.method == 3) + out.pca <- horizFPCA(out, M, showplot=F) + + cumm_coef = cumsum(out.pca$latent)/sum(out.pca$latent) + no = which(cumm_coef >= pc) + no = no[1] + + lam = 1/out.pca$latent[1:no] + Sigma = diag(lam) + eta = out.pca$coef[, 1:no] + eta_bar = apply(eta,2,sum) + + # compute test statistic + Sn <- rep(0, N1) + for (j in (2:N1)) { + tmp_eta = eta[1:j,] + tmp = apply(tmp_eta,2,sum)-dim(tmp_eta)[1]*eta_bar + + Sn[j] = 1/N1 * t(tmp) %*% Sigma %*% tmp + } + k.star <- min(which(Sn == max(Sn))) + Tn <- mean(Sn) + + # compute distribution + Values = rep(0, d) + for (i in (1:d)){ + B_tmp = matrix(0,no, N1) + for (j in (1:no)){ + B_tmp[j,] = BBridge(0, 0, 0, 1, N1-1)^2 } - - out <- time_warping(f, time, parallel = T, showplot=F) - - # Calculate PCA ----------------------------------------------------------- - if (pca.method == 1) - out.pca <- jointFPCA(out, N1, showplot=F, C=NULL) - if (pca.method == 2) - out.pca <- vertFPCA(out, N1, showplot=F) - if (pca.method == 3) - out.pca <- horizFPCA(out, N1, showplot=F) - - cumm_coef = cumsum(out.pca$latent)/sum(out.pca$latent) - no = which(cumm_coef >= pc) - no = no[1] - - lam = 1/out.pca.latent[1:no] - Sigma = diag(lam) - eta = out.pca$coef[, 1:no] - eta_bar = apply(eta,2,sum) - - # compute test statisitc - Sn <- rep(0, N1) - for (j in (2:N)) { - tmp_eta = eta[1:i,] - tmp = apply(tmp_eta,2,sum)-dim(tmp_eta)[1]*eta_bar - - Sn[j] = 1/N1 * t(tmp) %*% Sigma %*% tmp + Values[i] = mean(apply(B_tmp,2,sum)) + } + + z <- Tn <= Values + p <- length(z[z == TRUE]) / length(z) + dat.a <- f[, 1:k.star] + dat.b <- f[, (k.star + 1):N1] + warp.a <- out$warping_functions[,1:k.star] + warp.b <- out$warping_functions[,(k.star + 1):N1] + mean.a <- rowMeans(out$fn[,1:k.star]) + mean.b <- rowMeans(out$fn[,(k.star + 1):N1]) + out.gam.mean <- SqrtMean(warp.a) + warp_mean_a = out.gam.mean$gam_mu + out.gam.mean <- SqrtMean(warp.b) + warp_mean_b = out.gam.mean$gam_mu + delta <- mean.a - mean.b + + # Plot + delta <- mean.a - mean.b + if (showplot == TRUE) { + par(mfrow = c(1, 3)) + matplot(f, type = "l", col = "grey", main = "Functional Data", ylab = "values") + for (i in 1:ncol(dat.b)) { + lines(dat.b[, i], col = "pink") } - k.star <- min(which(Sn == max(Sn))) - Tn <- mean(Sn) - - # compute distribution - Values = rep(0, d) - for (i in (1:d)){ - B_tmp = matrix(0,no, N1) - for (j in (1:no)){ - B_tmp[j,] = BBridge(0, 0, 0, 1, N1-1)^2 - } - Values[i] = mean(apply(B_tmp,2,sum)) + for (i in 1:ncol(dat.a)) { + lines(dat.a[, i], col = "lightblue") } + lines(mean.a, col = "blue") + lines(mean.b, col = "red") + legend("topleft", c("before", "after"), col = c("blue", "red"), lty = c(1, 1), cex = 0.5) - z <- Tn <= Values - p <- length(z[z == TRUE]) / length(z) - dat.a <- f[, 1:k.star] - dat.b <- f[, (k.star + 1):N1] - warp.a <- out$warping_functions[,1:k.star] - warp.b <- out$warping_functions[,(k.star + 1):N1] - mean.a <- rowMeans(out$fn[,1:k.star]) - mean.b <- rowMeans(out$fn[,(k.star + 1):N1]) - out.gam.mean <- SqrtMean(warp.a) - warp_mean_a = out.gam.mean$gam_mu - out.gam.mean <- SqrtMean(warp.b) - warp_mean_b = out.gam.mean$gam_mu - delta <- mean.a - mean.b - - # Plot - delta <- mean.a - mean.b - if (showplot == TRUE) { - par(mfrow = c(1, 3)) - matplot(f, type = "l", col = "grey", main = "Functional Data", ylab = "values") - for (i in 1:ncol(dat.b)) { - lines(dat.b[, i], col = "pink") - } - for (i in 1:ncol(dat.a)) { - lines(dat.a[, i], col = "lightblue") - } - lines(mean.a, col = "blue") - lines(mean.b, col = "red") - legend("topleft", c("before", "after"), col = c("blue", "red"), lty = c(1, 1), cex = 0.5) - - plot(delta, type = "l", main = "Estimated Change Function", ylab = "values") - - matplot(out$warping_functions, type = "l", col = "grey", main = "Warping Functions", ylab = "values") - for (i in 1:ncol(warp.b)) { - lines(warp.b[, i], col = "pink") - } - for (i in 1:ncol(warp.a)) { - lines(warp.a[, i], col = "lightblue") - } - lines(warp_mean_a, col = "blue") - lines(warp_mean_b, col = "red") - legend("topleft", c("before", "after"), col = c("blue", "red"), lty = c(1, 1), cex = 0.5) - } + plot(delta, type = "l", main = "Estimated Change Function", ylab = "values") - out <- list( - pvalue = p, change = k.star, - DataBefore = dat.a, DataAfter = dat.b, - MeanBefore = mean.a, MeanAfter = mean.b, - WarpingMeanBefore = warp_mean_a, WarpingMeanAfter = warp_mean_b, - WarpingBefore = warp.a, WarpingAfter = warp.b, - change_fun = delta, Sn = Sn - ) - - return(out) + matplot(out$warping_functions, type = "l", col = "grey", main = "Warping Functions", ylab = "values") + for (i in 1:ncol(warp.b)) { + lines(warp.b[, i], col = "pink") + } + for (i in 1:ncol(warp.a)) { + lines(warp.a[, i], col = "lightblue") + } + lines(warp_mean_a, col = "blue") + lines(warp_mean_b, col = "red") + legend("topleft", c("before", "after"), col = c("blue", "red"), lty = c(1, 1), cex = 0.5) + } + + out <- list( + pvalue = p, change = k.star, + DataBefore = dat.a, DataAfter = dat.b, + MeanBefore = mean.a, MeanAfter = mean.b, + WarpingMeanBefore = warp_mean_a, WarpingMeanAfter = warp_mean_b, + WarpingBefore = warp.a, WarpingAfter = warp.b, + change_fun = delta, Sn = Sn + ) + + return(out) } diff --git a/man/elastic_amp_change_ff.Rd b/man/elastic_amp_change_ff.Rd index 9689080..63c5af3 100644 --- a/man/elastic_amp_change_ff.Rd +++ b/man/elastic_amp_change_ff.Rd @@ -42,7 +42,7 @@ Returns a list object containing \item{WarpingMeanBefore}{mean warping function before changepoint} \item{WarpingMeanAfter}{mean warping function after changepoint} \item{change_fun}{amplitude change function} -\item{Sn}{test statisitc values} +\item{Sn}{test statistic values} \item{mu}{mean srsfs} \item{mu_f}{mean functions} } diff --git a/man/elastic_change_fpca.Rd b/man/elastic_change_fpca.Rd index fa17406..817431d 100644 --- a/man/elastic_change_fpca.Rd +++ b/man/elastic_change_fpca.Rd @@ -49,7 +49,7 @@ Returns a list object containing \item{WarpingMeanBefore}{mean warping function before changepoint} \item{WarpingMeanAfter}{mean warping function after changepoint} \item{change_fun}{amplitude change function} -\item{Sn}{test statisitc values} +\item{Sn}{test statistic values} } \description{ This function identifies changepoints using a functional PCA diff --git a/man/elastic_ph_change_ff.Rd b/man/elastic_ph_change_ff.Rd index 0b1fa2e..9deefb5 100644 --- a/man/elastic_ph_change_ff.Rd +++ b/man/elastic_ph_change_ff.Rd @@ -42,7 +42,7 @@ Returns a list object containing \item{WarpingMeanBefore}{mean warping function before changepoint} \item{WarpingMeanAfter}{mean warping function after changepoint} \item{change_fun}{amplitude change function} -\item{Sn}{test statisitc values} +\item{Sn}{test statistic values} \item{mu}{mean shooting vectors} } \description{