# R/mfssa.R In Rfssa: Functional Singular Spectrum Analysis

#### Defines functions mfssa

```# Embedding and decomposition stages of multivariate functional singular spectrum analysis.
mfssa <- function(Y, L, ntriples) {
# get c plus plus code
p <- length(Y@C)
Y_d <- matrix(data = 0, nrow = 1, ncol = p)
for (j in 1:p) Y_d[1, j] <- ncol(Y@B[[j]])
d <- L * matrix(c(0, Y_d), nrow = 1L, ncol = (p + 1L))
N <- ncol(Y@C[[1]])
C_tilde <- list()
G_1 <- list()
# get inner product matrices
for (i in 1:p) {
if (ncol(Y@grid[[i]]) == 1) {
C_tilde[[i]] <- t(onedG(A = Y@B[[i]] %*% Y@C[[i]], B = Y@B[[i]], grid = Y@grid[[i]])) # old B
G_1[[i]] <- t(onedG(A = Y@B[[i]], B = Y@B[[i]], grid = Y@grid[[i]])) # old A
} else {
C_tilde[[i]] <- t(twodG(A = Y@B[[i]] %*% Y@C[[i]], B = Y@B[[i]], grid = Y@grid[[i]]))
G_1[[i]] <- t(twodG(A = Y@B[[i]], B = Y@B[[i]], grid = Y@grid[[i]]))
}
}
# Find the proper inner product matrices for j_k variables
d_tilde <- sum(d) / L
K <- N - L + 1L
shifter <- matrix(nrow = 2, ncol = (p + 1L), data = 0L)
shifter[, 2L] <- c(1L, d[2L])
if (p > 1L) {
for (i in 2L:p) {
shifter[1L, i + 1L] <- shifter[2L, i] + 1L
shifter[2L, i + 1L] <- shifter[2L, i] + d[i + 1L]
}
}
# Calculating Variance/Covariance Structure
S0 <- SSM(K, L, d_tilde, p, C_tilde, shifter)
# Calculating Gram Matrix
H <- CalculateInverse(Gramm(K, L, p, d_tilde, G_1, shifter, d))
# Calculating Eigen Triples
Q <- eigs(AtimesB(H, S0), ntriples)
# Returning results
Q\$values <- Re(Q\$values)
Q\$vectors <- Re(Q\$vectors)
coefs0 <- Q\$vectors
p_c <- list()
values <- Q\$values[1L:ntriples]
out <- list()
for (i in 1L:(ntriples)) {
my_pcs <- list(NA)
for (j in 1L:p) {
my_pcs[[j]] <- Y@B[[j]] %*% Cofmat((d[j + 1L] / L), L, coefs0[(shifter[1L, (j + 1L)]:shifter[2L, (j + 1L)]), i])
}
out[[i]] <- my_pcs
}
out\$values <- values
out\$L <- L
out\$N <- N
out\$Y <- Y
out\$RVectrs <- mV(out, ntriples)
return(out)
}
```

## Try the Rfssa package in your browser

Any scripts or data that you put into this service are public.

Rfssa documentation built on Jan. 10, 2022, 1:07 a.m.