popsimA: Simulation time series data for population A

In BrainCon: Inference the Partial Correlations Based on Time Series Data

Simulation time series data for population A

Description

A dataset containing values of 10 interested variables of 20 subjects over 50 periods.

Usage

popsimA

Format

An object of class array of dimension 50 x 10 x 20.

See Also

popsimB.

Examples

## Generated by the following R codes
set.seed(1234)
n = 50; p = 10; m1 = 20; m2 = 10
Precision1 = Precision2 = diag(rep(1, p))    # generate Precision matrix for population
for (i in 1 : (p - 1)){
  temp1 = ifelse(i > 2 * p / 3, -0.2, 0.4)
  temp2 = ifelse(i < p / 3, 0.4, -0.2)
  Precision1[i, i + 1] = Precision1[i + 1, i] = temp1
  Precision2[i, i + 1] = Precision2[i + 1, i] = temp2
}
# R1=-cov2cor(Precision1) + diag(rep(2, p))  # real partial correlation matrix
# R2=-cov2cor(Precision2) + diag(rep(2, p))
Index = matrix(0, p, p)                      # generate covariance matrix for each subject
for (i in 1 : p){
  for (j in 1 : p){
    if (i != j & abs(i - j) <= 3) Index[i, j] = 1
  }
}
SigmaAll1 = array(dim = c(p, p, m1))
SigmaAll2 = array(dim = c(p, p, m2))
for (sub in 1 : m1){
  RE = matrix(rnorm(p^2, 0, sqrt(2) * 0.05), p, p) * Index
  RE1 = (RE + t(RE)) / 2
  PrecisionInd = Precision1 + RE1
  SigmaAll1[, , sub] = solve(PrecisionInd)
}
for (sub in 1 : m2){
  RE = matrix(rnorm(p^2, 0, sqrt(2) * 0.15), p, p) * Index
  RE1 = (RE + t(RE)) / 2
  PrecisionInd = Precision2 + RE1
  SigmaAll2[, , sub] = solve(PrecisionInd)
}
rho = 0.3                                    # generate observed time series data
y1 = array(dim = c(n, p, m1))
y2 = array(dim = c(n, p, m2))
for (sub in 1 : m1){
  SigmaInd1 = SigmaAll1[, , sub]
  ytemp = matrix(0, n, p)
  Epsilon = MASS::mvrnorm(n, rep(0, p), SigmaInd1)
  ytemp[1, ] = Epsilon[1, ]
  for (i in 2 : n){
    ytemp[i, ] = rho * ytemp[i - 1, ] + sqrt(1 - rho^2) * Epsilon[i, ]
  }
  y1[, , sub] = ytemp
}
for (sub in 1 : m2){
  SigmaInd2 = SigmaAll2[, , sub]
  Xtemp = matrix(0, n, p)
  Epsilon = MASS::mvrnorm(n, rep(0, p), SigmaInd2)
  ytemp[1, ] = Epsilon[1, ]
  for (i in 2 : n){
    ytemp[i, ] = rho * ytemp[i - 1, ] + sqrt(1 - rho^2) * Epsilon[i, ]
  }
  y2[, , sub] = ytemp
}
popsimA = y1
popsimB = y2

BrainCon documentation built on April 21, 2022, 5:08 p.m.