NetFIntegral {frechet}

Generalized Fréchet integrals of network

Description

Calculating generalized Fréchet integrals of networks (equipped with Frobenius norm of adjacency matrices with zero diagonal elements and non negative off diagonal elements.)

Usage

NetFIntegral(phi, t_out, X, U)

Arguments

phi

An eigenfunction along which we want to project the network

t_out

Support of phi

X

A three dimensional array of dimension length(t_out) x m x m, where X[i,,] is an m x m network adjacency matrix. The diagonal elements of adjacency matrices are zero and the off diagonal entries lie between zero and U.

U

Upper bound of off-diagonal entries

Value

A list of the following:

f

An adjacency matrix which corresponds to the Fréchet integral of X along phi

References

Dubey, P., & Müller, H. G. (2020). Functional models for time‐varying random objects. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 82(2), 275-327.

Examples

set.seed(5)
n <- 100
N <- 50
t_out <- seq(0,1,length.out = N)
library(mpoly)
p2 <- as.function(mpoly::jacobi(2,4,3),silent=TRUE)
p4 <- as.function(mpoly::jacobi(4,4,3),silent=TRUE)
p6 <- as.function(mpoly::jacobi(6,4,3),silent=TRUE)

# first three eigenfunctions
phi1 <- function(t){
p2(2*t-1)*t^(1.5)*(1-t)^2 / (integrate(function(x) p2(2*x-1)^2*x^(3)*(1-x)^4,0,1))$value^(1/2)
}
phi2 <- function(t){
p4(2*t-1)*t^(1.5)*(1-t)^2 / (integrate(function(x) p4(2*x-1)^2*x^(3)*(1-x)^4,0,1))$value^(1/2)
}
phi3 <- function(t){
p6(2*t-1)*t^(1.5)*(1-t)^2 / (integrate(function(x) p6(2*x-1)^2*x^(3)*(1-x)^4,0,1))$value^(1/2)
}

# random component of adjacency matrices
P12 <- 0.1 ## edge between compunities
Score <- matrix(runif(n*4), nrow = n)
# edge within first community
P1_vec <- 0.5 + 0.4*Score[,1] %*% t(phi1(t_out)) + 0.1*Score[,2] %*% t(phi3(t_out)) 
# edge within second community
P2_vec <- 0.5 + 0.3*Score[,3] %*% t(phi2(t_out)) + 0.1*Score[,4] %*% t(phi3(t_out)) 

# create Network edge matrix
N_net1 <- 5 # first community number
N_net2 <- 5 # second community number

# I: four dimension array of n x n matrix of squared distances between the time point u 
# of the ith process and process and the time point v of the jth object process,
# e.g.: I[i,j,u,v] <- d_F^2(X_i(u) X_j(v)).
I <- array(0, dim = c(n,n,N,N))
for(u in 1:N){
  for(v in 1:N){
   #frobenius norm between two adjcent matrix
    I[,,u,v] <- outer(P1_vec[,u], P1_vec[,v], function(a1, a2) (a1-a2)^2*(N_net1^2-N_net1)) +
      outer(P2_vec[,u], P2_vec[,v], function(a1, a2) (a1-a2)^2*(N_net2^2-N_net2))
  }
}


# check ObjCov work
Cov_result <- ObjCov(t_out, I, 3, smooth=FALSE)
Cov_result$lambda  # 0.266 0.15 0.04

# sum((Cov_result$phi[,1] - phi1(t_out))^2) / sum(phi1(t_out)^2)
# sum((Cov_result$phi[,2] - phi2(t_out))^2) / sum(phi2(t_out)^2)
# sum((Cov_result$phi[,3] - phi3(t_out))^2) / sum(phi3(t_out)^2)

# e.g. subj 2
subj <- 2
# X_mat is the network for varying times with X[i,,] is the adjacency matrices 
# for the ith time point
X_mat <- array(0, c(N,(N_net1+N_net2), (N_net1+N_net2)))
for(i in 1:N){
  # edge between communities is P12
  Mat <- matrix(P12, nrow = (N_net1+N_net2), ncol = (N_net1+N_net2)) 
  # edge within the first communitiy is P1
  Mat[1:N_net1, 1:N_net1] <- P1_vec[subj, i] 
  # edge within the second community is P2
  Mat[(N_net1+1):(N_net1+N_net2), (N_net1+1):(N_net1+N_net2)] <- P2_vec[subj, i] 
  diag(Mat) <- 0 #diagonal element is 0
  X_mat[i,,] <- Mat
}
# output the functional principal network(adjacency matrice) of the second eigenfunction
NetFIntegral(Cov_result$phi[,2], t_out, X_mat, 2)