summary_sys: Summary of Correlated Systems of Variables
In SimRepeat: Simulation of Correlated Systems of Equations with Multiple Variable Types
Description
Usage
Arguments
Value
References
See Also
Examples
Description
This function summarizes the results of nonnormsys, corrsys, or
corrsys2. The inputs are either the simulated variables or inputs for those functions. See their
documentation for more information. If only selected descriptions are desired, keep the non-relevant parameter inputs at their
defaults. For example, if only a description of the error terms are desired, error_type = "non_mix", and
method = "Polynomial", specify E, M, method, means, vars, skews, skurts, fifths, sixths, corr.e.
Usage
1
2
3
4
5
6
7
8
9
10
11
summary_sys(Y = NULL, E = NULL, E_mix = NULL, X = list(),
X_all = list(), M = NULL, method = c("Fleishman", "Polynomial"),
means = list(), vars = list(), skews = list(), skurts = list(),
fifths = list(), sixths = list(), mix_pis = list(), mix_mus = list(),
mix_sigmas = list(), mix_skews = list(), mix_skurts = list(),
mix_fifths = list(), mix_sixths = list(), marginal = list(),
support = list(), lam = list(), p_zip = list(), size = list(),
prob = list(), mu = list(), p_zinb = list(), corr.x = list(),
corr.e = NULL, U = list(), U_all = list(), rand.int = c("none",
"non_mix", "mix"), rand.tsl = c("none", "non_mix", "mix"),
corr.u = list(), rmeans2 = list(), rvars2 = list())
Arguments
Y
the matrix of outcomes simulated with corrsys or corrsys2
E
the matrix of continuous non-mixture or components of mixture error terms
E_mix
the matrix of continuous mixture error terms
X
a list of length M where X[[p]] = cbind(X_cat(pj), X_cont(pj), X_comp(pj), X_pois(pj), X_nb(pj)); keep
X[[p]] = NULL if Y_p has no independent variables
X_all
a list of length M where X_all[[p]] contains all independent variables, interactions, and time for Y_p; keep
X_all[[p]] = NULL if Y_p has no independent variables
M
the number of dependent variables Y (outcomes); equivalently, the number of equations in the system
method
the PMT method used to generate all continuous variables, including independent variables (covariates), error terms, and random effects;
"Fleishman" uses Fleishman's third-order polynomial transformation and "Polynomial" uses Headrick's fifth-order transformation
means
if no random effects, a list of length M where means[[p]] contains a vector of means for the continuous independent variables
in equation p with non-mixture (X_{cont}) or mixture (X_{mix}) distributions and for the error terms (E);
order in vector is X_{cont}, X_{mix}, E
if there are random effects, a list of length M + 1 if the effects are the same across equations or 2 * M if they differ;
where means[M + 1] or means[(M + 1):(2 * M)] are vectors of means for all random effects with continuous non-mixture or mixture distributions;
order in vector is 1st random intercept U_0 (if rand.int != "none"), 2nd random time slope U_1 (if rand.tsl != "none"),
3rd other random slopes with non-mixture distributions U_{cont}, 4th other random slopes with mixture distributions U_{mix}
vars
a list of same length and order as means containing vectors of variances for the continuous variables, error terms, and any random effects
skews
if no random effects, a list of length M where skews[[p]] contains a vector of skew values for the continuous independent variables
in equation p with non-mixture (X_{cont}) distributions and for E if error_type = "non_mix"; order in vector is X_{cont}, E
if there are random effects, a list of length M + 1 if the effects are the same across equations or 2 * M if they differ;
where skews[M + 1] or skews[(M + 1):(2 * M)] are vectors of skew values for all random effects with continuous non-mixture distributions;
order in vector is 1st random intercept U_0 (if rand.int = "non_mix"), 2nd random time slope U_1 (if rand.tsl = "non_mix"),
3rd other random slopes with non-mixture distributions U_{cont}
skurts
a list of same length and order as skews containing vectors of standardized kurtoses (kurtosis - 3) for the continuous variables,
error terms, and any random effects with non-mixture distributions
fifths
a list of same length and order as skews containing vectors of standardized fifth cumulants for the continuous variables,
error terms, and any random effects with non-mixture distributions; not necessary for method = "Fleishman"
sixths
a list of same length and order as skews containing vectors of standardized sixth cumulants for the continuous variables,
error terms, and any random effects with non-mixture distributions; not necessary for method = "Fleishman"
mix_pis
list of length M, M + 1 or 2 * M, where mix_pis[1:M] are for X_{cont}, E (if error_type = "mix") and
mix_pis[M + 1] or mix_pis[(M + 1):(2 * M)] are for mixture U; use mix_pis[[p]] = NULL if equation p has no continuous mixture terms
if error_type = "non_mix" and there are only random effects (i.e., length(corr.x) = 0), use mix_pis[1:M] = NULL so that
mix_pis[M + 1] or mix_pis[(M + 1):(2 * M)] describes the mixture U;
mix_pis[[p]][[j]] is a vector of mixing probabilities of the component distributions for X_{mix(pj)}, the j-th mixture covariate for outcome Y_p;
the last vector in mix_pis[[p]] is for E_p (if error_type = "mix"); components should be ordered as in corr.x
mix_pis[[M + p]][[j]] is a vector of mixing probabilities of the component distributions for U_{(pj)}, the j-th random effect with a mixture distribution
for outcome Y_p; order is 1st random intercept (if rand.int = "mix"), 2nd random time slope (if rand.tsl = "mix"),
3rd other random slopes with mixture distributions; components should be ordered as in corr.u
mix_mus
list of same length and order as mix_pis;
mix_mus[[p]][[j]] is a vector of means of the component distributions for X_{mix(pj)},
the last vector in mix_mus[[p]] is for E_p (if error_type = "mix")
mix_mus[[p]][[j]] is a vector of means of the component distributions for U_{mix(pj)}
mix_sigmas
list of same length and order as mix_pis;
mix_sigmas[[p]][[j]] is a vector of standard deviations of the component distributions for X_{mix(pj)},
the last vector in mix_sigmas[[p]] is for E_p (if error_type = "mix")
mix_sigmas[[p]][[j]] is a vector of standard deviations of the component distributions for U_{mix(pj)}
mix_skews
list of same length and order as mix_pis;
mix_skews[[p]][[j]] is a vector of skew values of the component distributions for X_{mix(pj)},
the last vector in mix_skews[[p]] is for E_p (if error_type = "mix")
mix_skews[[p]][[j]] is a vector of skew values of the component distributions for U_{mix(pj)}
mix_skurts
list of same length and order as mix_pis;
mix_skurts[[p]][[j]] is a vector of standardized kurtoses of the component distributions for X_{mix(pj)},
the last vector in mix_skurts[[p]] is for E_p (if error_type = "mix")
mix_skurts[[p]][[j]] is a vector of standardized kurtoses of the component distributions for U_{mix(pj)}
mix_fifths
list of same length and order as mix_pis; not necessary for method = "Fleishman";
mix_fifths[[p]][[j]] is a vector of standardized fifth cumulants of the component distributions for X_{mix(pj)},
the last vector in mix_fifths[[p]] is for E_p (if error_type = "mix")
mix_fifths[[p]][[j]] is a vector of standardized fifth cumulants of the component distributions for U_{mix(pj)}
mix_sixths
list of same length and order as mix_pis; not necessary for method = "Fleishman";
mix_sixths[[p]][[j]] is a vector of standardized sixth cumulants of the component distributions for X_{mix(pj)},
the last vector in mix_sixths[[p]] is for E_p (if error_type = "mix")
mix_sixths[[p]][[j]] is a vector of standardized sixth cumulants of the component distributions for U_{mix(pj)}
marginal
a list of length M, with the p-th component a list of cumulative probabilities for the ordinal variables
associated with outcome Y_p (use marginal[[p]] = NULL if outcome Y_p has no ordinal variables);
marginal[[p]][[j]] is a vector of the cumulative probabilities defining the marginal distribution of X_{ord(pj)},
the j-th ordinal variable for outcome Y_p; if the variable can take r values, the vector will contain r - 1 probabilities
(the r-th is assumed to be 1); for binary variables, the probability is the probability of the 1st category, which has the smaller support value;
length(marginal[[p]]) can differ across outcomes; the order should be the same as in corr.x
support
a list of length M, with the p-th component a list of support values for the ordinal variables associated
with outcome Y_p; use support[[p]] = NULL if outcome Y_p has no ordinal variables;
support[[p]][[j]] is a vector of the support values defining the marginal distribution of X_{ord(pj)},
the j-th ordinal variable for outcome Y_p; if not provided, the default for r categories is 1, ..., r
lam
list of length M, p-th component a vector of lambda (means > 0) values for Poisson variables for outcome Y_p
(see stats::dpois); order is 1st regular Poisson and 2nd zero-inflated Poisson; use lam[[p]] = NULL if outcome Y_p has no Poisson variables;
length(lam[[p]]) can differ across outcomes; the order should be the same as in corr.x
p_zip
a list of vectors of probabilities of structural zeros (not including zeros from the Poisson distribution) for the
zero-inflated Poisson variables (see VGAM::dzipois); if p_zip = 0, Y_{pois} has a regular Poisson
distribution; if p_zip is in (0, 1), Y_{pois} has a zero-inflated Poisson distribution;
if p_zip is in (-(exp(lam) - 1)^(-1), 0), Y_{pois} has a zero-deflated Poisson distribution and p_zip
is not a probability; if p_zip = -(exp(lam) - 1)^(-1), Y_{pois} has a positive-Poisson distribution
(see VGAM::dpospois); order is 1st regular Poisson and 2nd zero-inflated Poisson;
if a single number, all Poisson variables given this value; if a vector of length M, all Poisson variables in equation p
given p_zip[p]; otherwise, missing values are set to 0 and ordered 1st
size
list of length M, p-th component a vector of size parameters for the Negative Binomial variables for outcome Y_p
(see stats::nbinom); order is 1st regular NB and 2nd zero-inflated NB; use size[[p]] = NULL if outcome Y_p has no Negative Binomial variables;
length(size[[p]]) can differ across outcomes; the order should be the same as in corr.x
prob
list of length M, p-th component a vector of success probabilities for the Negative Binomial variables for outcome Y_p
(see stats::nbinom); order is 1st regular NB and 2nd zero-inflated NB; use prob[[p]] = NULL if outcome Y_p has no Negative Binomial variables;
length(prob[[p]]) can differ across outcomes; the order should be the same as in corr.x
mu
list of length M, p-th component a vector of mean values for the Negative Binomial variables for outcome Y_p
(see stats::nbinom); order is 1st regular NB and 2nd zero-inflated NB; use mu[[p]] = NULL if outcome Y_p has no Negative Binomial variables;
length(mu[[p]]) can differ across outcomes; the order should be the same as in corr.x; for zero-inflated NB variables,
this refers to the mean of the NB distribution (see VGAM::dzinegbin)
(*Note: either prob or mu should be supplied for all Negative Binomial variables, not a mixture)
p_zinb
a vector of probabilities of structural zeros (not including zeros from the NB distribution) for the zero-inflated NB variables
(see VGAM::dzinegbin); if p_zinb = 0, Y_{nb} has a regular NB distribution;
if p_zinb is in (-prob^size/(1 - prob^size), 0), Y_{nb} has a zero-deflated NB distribution and p_zinb
is not a probability; if p_zinb = -prob^size/(1 - prob^size), Y_{nb} has a positive-NB distribution (see
VGAM::dposnegbin); order is 1st regular NB and 2nd zero-inflated NB;
if a single number, all NB variables given this value; if a vector of length M, all NB variables in equation p
given p_zinb[p]; otherwise, missing values are set to 0 and ordered 1st
corr.x
list of length M, each component a list of length M; corr.x[[p]][[q]] is matrix of correlations
for independent variables in equations p (X_{(pj)} for outcome Y_p) and q (X_{(qj)} for outcome Y_q);
order: 1st ordinal (same order as in marginal), 2nd continuous non-mixture (same order as in skews),
3rd components of continuous mixture (same order as in mix_pis), 4th regular Poisson, 5th zero-inflated Poisson (same order as in lam),
6th regular NB, and 7th zero-inflated NB (same order as in size);
if p = q, corr.x[[p]][[q]] is a correlation matrix with nrow(corr.x[[p]][[q]]) = # X_{(pj)} for outcome Y_p;
if p != q, corr.x[[p]][[q]] is a non-symmetric matrix of correlations where rows correspond to covariates for Y_p
so that nrow(corr.x[[p]][[q]]) = # X_{(pj)} for outcome Y_p and
columns correspond to covariates for Y_q so that ncol(corr.x[[p]][[q]]) = # X_{(qj)} for outcome Y_q;
use corr.x[[p]][[q]] = NULL if equation q has no X_{(qj)}; use corr.x[[p]] = NULL if equation p has no X_{(pj)}
corr.e
correlation matrix for continuous non-mixture or components of mixture error terms
U
a list of length M of continuous non-mixture and components of mixture random effects
U_all
a list of length M of continuous non-mixture and mixture random effects
rand.int
"none" (default) if no random intercept term for all outcomes, "non_mix" if all random intercepts have a continuous
non-mixture distribution, "mix" if all random intercepts have a continuous mixture distribution;
also can be a vector of length M containing a combination (i.e., c("non_mix", "mix", "none") if the 1st has a non-mixture
distribution, the 2nd has a mixture distribution, and 3rd outcome has no random intercept)
rand.tsl
"none" (default) if no random slope for time for all outcomes, "non_mix" if all random time slopes have a
continuous non-mixture distribution, "mix" if all random time slopes have a continuous mixture distribution; also can
be a vector of length M as in rand.int
corr.u
if the random effects are the same variables across equations, a matrix of correlations for U;
if the random effects are different variables across equations, a list of length M, each component a list of length M;
corr.u[[p]][[q]] is matrix of correlations for random effects in equations p (U_{(pj)} for outcome Y_p) and
q (U_{(qj)} for outcome Y_q);
if p = q, corr.u[[p]][[q]] is a correlation matrix with nrow(corr.u[[p]][[q]]) = # U_{(pj)} for outcome Y_p;
if p != q, corr.u[[p]][[q]] is a non-symmetric matrix of correlations where rows correspond to U_{(pj)} for Y_p
so that nrow(corr.u[[p]][[q]]) = # U_{(pj)} for outcome Y_p and
columns correspond to U_{(qj)} for Y_q so that ncol(corr.u[[p]][[q]]) = # U_{(qj)} for outcome Y_q;
the number of random effects for Y_p is taken from nrow(corr.u[[p]][[1]]) so that if there should be random effects,
there must be entries for corr.u;
use corr.u[[p]][[q]] = NULL if equation q has no U_{(qj)}; use corr.u[[p]] = NULL if equation p has no U_{(pj)};
correlations are specified in terms of components of mixture variables (if present);
order is 1st random intercept (if rand.int != "none"), 2nd random time slope (if rand.tsl != "none"),
3rd other random slopes with non-mixture distributions, 4th other random slopes with mixture distributions
rmeans2
a list returned from corrsys or corrsys2 which has the non-mixture and component means ordered according to
types of random intercept and time slope
rvars2
a list returned like rmeans
Value
A list with the following components:
cont_sum_y a data.frame summarizing the simulated distributions of the Y_p,
cont_sum_e a data.frame summarizing the simulated distributions of the non-mixture or components of mixture E_p,
target_sum_e a data.frame summarizing the target distributions of the non-mixture or components of mixture E_p,
mix_sum_e a data.frame summarizing the simulated distributions of the mixture E_p,
target_mix_e a data.frame summarizing the target distributions of the mixture E_p,
rho.y correlation matrix of dimension M x M for Y_p
rho.e correlation matrix for the non-mixture or components of mixture E_p
rho.emix correlation matrix for the mixture E_p
rho.ye matrix with correlations between Y_p (rows) and the non-mixture or components of mixture E_p (columns)
rho.yemix matrix with correlations between Y_p (rows) and the mixture E_p (columns)
sum_xall a data.frame summarizing X_all without the Time variable,
rho.yx a list of length M, where rho.yx[[p]] is matrix of correlations
between Y (rows) and X[[p]] = {X_ord(pj), X_cont(pj), X_comp(pj), X_pois(pj), X_nb(pj)} (columns)
rho.yxall a list of length M, where rho.yx[[p]] is matrix of correlations
between Y (rows) and X_all[[p]] (columns) not including Time
rho.x a list of length M of lists of length M where
rho.x[[p]][[q]] = cor(cbind(X[[p]], X[[q]])) if p!= q or
rho.x[[p]][[q]] = cor(X[[p]])) if p = q, where X[[p]] = {X_ord(pj), X_cont(pj), X_comp(pj), X_pois(pj), X_nb(pj)}
rho.xall a list of length M of lists of length M where
rho.xall[[p]][[q]] = cor(cbind(X_all[[p]], X_all[[q]])) if p!= q or
rho.xall[[p]][[q]] = cor(X_all[[p]])) if p = q, not including Time
maxerr a list of length M containing a vector of length M with the maximum correlation errors between outcomes,
maxerr[[p]]][q] = abs(max(corr.x[[p]][[q]] - rho.x[[p]][[q]]))
Additional components vary based on the type of simulated variables:
If ordinal variables are produced:
ord_sum_x a list where ord_sum_x[[j]] is a data.frame summarizing X_{ord(pj)} for all p = 1, ..., M
If continuous variables are produced:
cont_sum_x a data.frame summarizing the simulated distributions of the X_{cont(pj)} and X_comp(pj),
target_sum_x a data.frame summarizing the target distributions of the X_{cont(pj)} and X_comp(pj),
mix_sum_x a data.frame summarizing the simulated distributions of the X_{mix(pj)},
target_mix_x a data.frame summarizing the target distributions of the X_{mix(pj)}
If Poisson variables are produced:
pois_sum_x a data.frame summarizing the simulated distributions of the X_{pois(pj)}
If Negative Binomial variables are produced:
nb_sum_x a data.frame summarizing the simulated distributions of the X_{nb(pj)}
If random effects are produced:
cont_sum_u a data.frame summarizing the simulated distributions of the U_{cont(pj)} and U_{comp(pj)},
target_sum_u a data.frame summarizing the target distributions of the U_{cont(pj)} and U_{comp(pj)},
sum_uall a data.frame summarizing the simulated distributions of U_all,
mix_sum_u a data.frame summarizing the simulated distributions of the U_{mix(pj)},
target_mix_u a data.frame summarizing the target distributions of the U_{mix(pj)},
rho.u list of length M, each component a list of length M;
rho.u[[p]][[q]] = cor(cbind(U[[p]], U[[q]])) if p != q or rho.u[[p]][[q]] = cor(U[[p]])) if p = q
rho.uall list of length M, each component a list of length M;
rho.uall[[p]][[q]] = cor(cbind(U_all[[p]], U_all[[q]])) if p != q or rho.uall[[p]][[q]] = cor(U_all[[p]])) if p = q
maxerr_u list of length M containing a vector of length M with the maximum correlation errors for U between outcomes
maxerr_u[[p]]][q] = abs(max(corr.u[[p]][[q]] - rho.u[[p]][[q]]))
References
See references for SimRepeat.
See Also
Examples
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
M <- 3
B <- calc_theory("Beta", c(4, 1.5))
skews <- lapply(seq_len(M), function(x) B[3])
skurts <- lapply(seq_len(M), function(x) B[4])
fifths <- lapply(seq_len(M), function(x) B[5])
sixths <- lapply(seq_len(M), function(x) B[6])
Six <- lapply(seq_len(M), function(x) list(0.03))
corr.e <- matrix(c(1, 0.4, 0.4^2, 0.4, 1, 0.4, 0.4^2, 0.4, 1), M, M,
byrow = TRUE)
means <- lapply(seq_len(M), function(x) B[1])
vars <- lapply(seq_len(M), function(x) B[2]^2)
marginal <- list(0.3, 0.4, 0.5)
support <- lapply(seq_len(M), function(x) list(0:1))
corr.x <- list(list(matrix(1, 1, 1), matrix(0.4, 1, 1), matrix(0.4, 1, 1)),
list(matrix(0.4, 1, 1), matrix(1, 1, 1), matrix(0.4, 1, 1)),
list(matrix(0.4, 1, 1), matrix(0.4, 1, 1), matrix(1, 1, 1)))
betas <- list(0.5)
betas.t <- 1
betas.tint <- list(0.25)
Sys1 <- corrsys(10000, M, Time = 1:M, "Polynomial", "non_mix", means, vars,
skews, skurts, fifths, sixths, Six, marginal = marginal, support = support,
corr.x = corr.x, corr.e = corr.e, betas = betas, betas.t = betas.t,
betas.tint = betas.tint, quiet = TRUE)
Sum1 <- summary_sys(Sys1$Y, Sys1$E, E_mix = NULL, Sys1$X, Sys1$X_all, M,
"Polynomial", means, vars, skews, skurts, fifths, sixths,
marginal = marginal, support = support, corr.x = corr.x, corr.e = corr.e)
## Not run:
seed <- 276
n <- 10000
M <- 3
Time <- 1:M
# Error terms have a beta(4, 1.5) distribution with an AR(1, p = 0.4)
# correlation structure
B <- calc_theory("Beta", c(4, 1.5))
skews <- lapply(seq_len(M), function(x) B[3])
skurts <- lapply(seq_len(M), function(x) B[4])
fifths <- lapply(seq_len(M), function(x) B[5])
sixths <- lapply(seq_len(M), function(x) B[6])
Six <- lapply(seq_len(M), function(x) list(0.03))
error_type <- "non_mix"
corr.e <- matrix(c(1, 0.4, 0.4^2, 0.4, 1, 0.4, 0.4^2, 0.4, 1), M, M,
byrow = TRUE)
# 1 continuous mixture of Normal(-2, 1) and Normal(2, 1) for each Y
mix_pis <- lapply(seq_len(M), function(x) list(c(0.4, 0.6)))
mix_mus <- lapply(seq_len(M), function(x) list(c(-2, 2)))
mix_sigmas <- lapply(seq_len(M), function(x) list(c(1, 1)))
mix_skews <- lapply(seq_len(M), function(x) list(c(0, 0)))
mix_skurts <- lapply(seq_len(M), function(x) list(c(0, 0)))
mix_fifths <- lapply(seq_len(M), function(x) list(c(0, 0)))
mix_sixths <- lapply(seq_len(M), function(x) list(c(0, 0)))
mix_Six <- list()
Nstcum <- calc_mixmoments(mix_pis[[1]][[1]], mix_mus[[1]][[1]],
mix_sigmas[[1]][[1]], mix_skews[[1]][[1]], mix_skurts[[1]][[1]],
mix_fifths[[1]][[1]], mix_sixths[[1]][[1]])
means <- lapply(seq_len(M), function(x) c(Nstcum[1], B[1]))
vars <- lapply(seq_len(M), function(x) c(Nstcum[2]^2, B[2]^2))
# 1 binary variable for each Y
marginal <- lapply(seq_len(M), function(x) list(0.4))
support <- list(NULL, list(c(0, 1)), NULL)
# 1 Poisson variable for each Y
lam <- list(1, 5, 10)
# Y2 and Y3 are zero-inflated Poisson variables
p_zip <- list(NULL, 0.05, 0.1)
# 1 NB variable for each Y
size <- list(10, 15, 20)
prob <- list(0.3, 0.4, 0.5)
# either prob or mu is required (not both)
mu <- mapply(function(x, y) x * (1 - y)/y, size, prob, SIMPLIFY = FALSE)
# Y2 and Y3 are zero-inflated NB variables
p_zinb <- list(NULL, 0.05, 0.1)
# The 2nd (the normal mixture) variable is the same across Y
same.var <- 2
# Create the correlation matrix in terms of the components of the normal
# mixture
K <- 5
corr.x <- list()
corr.x[[1]] <- list(matrix(0.1, K, K), matrix(0.2, K, K), matrix(0.3, K, K))
diag(corr.x[[1]][[1]]) <- 1
# set correlation between components to 0
corr.x[[1]][[1]][2:3, 2:3] <- diag(2)
# set correlations with the same variable equal across outcomes
corr.x[[1]][[2]][, same.var] <- corr.x[[1]][[3]][, same.var] <-
corr.x[[1]][[1]][, same.var]
corr.x[[2]] <- list(t(corr.x[[1]][[2]]), matrix(0.35, K, K),
matrix(0.4, K, K))
diag(corr.x[[2]][[2]]) <- 1
corr.x[[2]][[2]][2:3, 2:3] <- diag(2)
corr.x[[2]][[2]][, same.var] <- corr.x[[2]][[3]][, same.var] <-
t(corr.x[[1]][[2]][same.var, ])
corr.x[[2]][[3]][same.var, ] <- corr.x[[1]][[3]][same.var, ]
corr.x[[2]][[2]][same.var, ] <- t(corr.x[[2]][[2]][, same.var])
corr.x[[3]] <- list(t(corr.x[[1]][[3]]), t(corr.x[[2]][[3]]),
matrix(0.5, K, K))
diag(corr.x[[3]][[3]]) <- 1
corr.x[[3]][[3]][2:3, 2:3] <- diag(2)
corr.x[[3]][[3]][, same.var] <- t(corr.x[[1]][[3]][same.var, ])
corr.x[[3]][[3]][same.var, ] <- t(corr.x[[3]][[3]][, same.var])
# The 2nd and 3rd variables of each Y are subject-level variables
subj.var <- matrix(c(1, 2, 1, 3, 2, 2, 2, 3, 3, 2, 3, 3), 6, 2, byrow = TRUE)
int.var <- tint.var <- NULL
betas.0 <- 0
betas <- list(seq(0.5, 0.5 + (K - 2) * 0.25, 0.25))
betas.subj <- list(seq(0.5, 0.5 + (K - 2) * 0.1, 0.1))
betas.int <- list()
betas.t <- 1
betas.tint <- list(c(0.25, 0.5))
method <- "Polynomial"
# Check parameter inputs
checkpar(M, method, error_type, means, vars, skews, skurts, fifths, sixths,
Six, mix_pis, mix_mus, mix_sigmas, mix_skews, mix_skurts, mix_fifths,
mix_sixths, mix_Six, marginal, support, lam, p_zip, pois_eps = list(),
size, prob, mu, p_zinb, nb_eps = list(), corr.x, corr.yx = list(),
corr.e, same.var, subj.var, int.var, tint.var, betas.0, betas,
betas.subj, betas.int, betas.t, betas.tint)
# Simulated system using correlation method 1
N <- corrsys(n, M, Time, method, error_type, means, vars, skews, skurts,
fifths, sixths, Six, mix_pis, mix_mus, mix_sigmas, mix_skews, mix_skurts,
mix_fifths, mix_sixths, mix_Six, marginal, support, lam, p_zip, size,
prob, mu, p_zinb, corr.x, corr.e, same.var, subj.var, int.var, tint.var,
betas.0, betas, betas.subj, betas.int, betas.t, betas.tint, seed = seed,
use.nearPD = FALSE)
# Summarize the results
S <- summary_sys(N$Y, N$E, E_mix = NULL, N$X, N$X_all, M, method, means,
vars, skews, skurts, fifths, sixths, mix_pis, mix_mus, mix_sigmas,
mix_skews, mix_skurts, mix_fifths, mix_sixths, marginal, support, lam,
p_zip, size, prob, mu, p_zinb, corr.x, corr.e)
S$sum_xall
S$maxerr
## End(Not run)
SimRepeat documentation built on May 2, 2019, 9:32 a.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Description Usage Arguments Value References See Also Examples
Description
This function summarizes the results of nonnormsys, corrsys, or
corrsys2. The inputs are either the simulated variables or inputs for those functions. See their
documentation for more information. If only selected descriptions are desired, keep the non-relevant parameter inputs at their
defaults. For example, if only a description of the error terms are desired, error_type = "non_mix", and
method = "Polynomial", specify E, M, method, means, vars, skews, skurts, fifths, sixths, corr.e.
Usage
1 2 3 4 5 6 7 8 9 10 11 | summary_sys(Y = NULL, E = NULL, E_mix = NULL, X = list(),
X_all = list(), M = NULL, method = c("Fleishman", "Polynomial"),
means = list(), vars = list(), skews = list(), skurts = list(),
fifths = list(), sixths = list(), mix_pis = list(), mix_mus = list(),
mix_sigmas = list(), mix_skews = list(), mix_skurts = list(),
mix_fifths = list(), mix_sixths = list(), marginal = list(),
support = list(), lam = list(), p_zip = list(), size = list(),
prob = list(), mu = list(), p_zinb = list(), corr.x = list(),
corr.e = NULL, U = list(), U_all = list(), rand.int = c("none",
"non_mix", "mix"), rand.tsl = c("none", "non_mix", "mix"),
corr.u = list(), rmeans2 = list(), rvars2 = list())
|
Arguments
Y |
the matrix of outcomes simulated with |
E |
the matrix of continuous non-mixture or components of mixture error terms |
E_mix |
the matrix of continuous mixture error terms |
X |
a list of length |
X_all |
a list of length |
M |
the number of dependent variables Y (outcomes); equivalently, the number of equations in the system |
method |
the PMT method used to generate all continuous variables, including independent variables (covariates), error terms, and random effects; "Fleishman" uses Fleishman's third-order polynomial transformation and "Polynomial" uses Headrick's fifth-order transformation |
means |
if no random effects, a list of length if there are random effects, a list of length |
vars |
a list of same length and order as |
skews |
if no random effects, a list of length if there are random effects, a list of length |
skurts |
a list of same length and order as |
fifths |
a list of same length and order as |
sixths |
a list of same length and order as |
mix_pis |
list of length
|
mix_mus |
list of same length and order as
|
mix_sigmas |
list of same length and order as
|
mix_skews |
list of same length and order as
|
mix_skurts |
list of same length and order as
|
mix_fifths |
list of same length and order as
|
mix_sixths |
list of same length and order as
|
marginal |
a list of length |
support |
a list of length |
lam |
list of length |
p_zip |
a list of vectors of probabilities of structural zeros (not including zeros from the Poisson distribution) for the
zero-inflated Poisson variables (see |
size |
list of length |
prob |
list of length |
mu |
list of length |
p_zinb |
a vector of probabilities of structural zeros (not including zeros from the NB distribution) for the zero-inflated NB variables
(see |
corr.x |
list of length |
corr.e |
correlation matrix for continuous non-mixture or components of mixture error terms |
U |
a list of length |
U_all |
a list of length |
rand.int |
"none" (default) if no random intercept term for all outcomes, "non_mix" if all random intercepts have a continuous
non-mixture distribution, "mix" if all random intercepts have a continuous mixture distribution;
also can be a vector of length |
rand.tsl |
"none" (default) if no random slope for time for all outcomes, "non_mix" if all random time slopes have a
continuous non-mixture distribution, "mix" if all random time slopes have a continuous mixture distribution; also can
be a vector of length |
corr.u |
if the random effects are the same variables across equations, a matrix of correlations for U;
if the random effects are different variables across equations, a list of length correlations are specified in terms of components of mixture variables (if present);
order is 1st random intercept (if |
rmeans2 |
a list returned from |
rvars2 |
a list returned like |
Value
A list with the following components:
cont_sum_y a data.frame summarizing the simulated distributions of the Y_p,
cont_sum_e a data.frame summarizing the simulated distributions of the non-mixture or components of mixture E_p,
target_sum_e a data.frame summarizing the target distributions of the non-mixture or components of mixture E_p,
mix_sum_e a data.frame summarizing the simulated distributions of the mixture E_p,
target_mix_e a data.frame summarizing the target distributions of the mixture E_p,
rho.y correlation matrix of dimension M x M for Y_p
rho.e correlation matrix for the non-mixture or components of mixture E_p
rho.emix correlation matrix for the mixture E_p
rho.ye matrix with correlations between Y_p (rows) and the non-mixture or components of mixture E_p (columns)
rho.yemix matrix with correlations between Y_p (rows) and the mixture E_p (columns)
sum_xall a data.frame summarizing X_all without the Time variable,
rho.yx a list of length M, where rho.yx[[p]] is matrix of correlations
between Y (rows) and X[[p]] = {X_ord(pj), X_cont(pj), X_comp(pj), X_pois(pj), X_nb(pj)} (columns)
rho.yxall a list of length M, where rho.yx[[p]] is matrix of correlations
between Y (rows) and X_all[[p]] (columns) not including Time
rho.x a list of length M of lists of length M where
rho.x[[p]][[q]] = cor(cbind(X[[p]], X[[q]])) if p!= q or
rho.x[[p]][[q]] = cor(X[[p]])) if p = q, where X[[p]] = {X_ord(pj), X_cont(pj), X_comp(pj), X_pois(pj), X_nb(pj)}
rho.xall a list of length M of lists of length M where
rho.xall[[p]][[q]] = cor(cbind(X_all[[p]], X_all[[q]])) if p!= q or
rho.xall[[p]][[q]] = cor(X_all[[p]])) if p = q, not including Time
maxerr a list of length M containing a vector of length M with the maximum correlation errors between outcomes,
maxerr[[p]]][q] = abs(max(corr.x[[p]][[q]] - rho.x[[p]][[q]]))
Additional components vary based on the type of simulated variables:
If ordinal variables are produced:
ord_sum_x a list where ord_sum_x[[j]] is a data.frame summarizing X_{ord(pj)} for all p = 1, ..., M
If continuous variables are produced:
cont_sum_x a data.frame summarizing the simulated distributions of the X_{cont(pj)} and X_comp(pj),
target_sum_x a data.frame summarizing the target distributions of the X_{cont(pj)} and X_comp(pj),
mix_sum_x a data.frame summarizing the simulated distributions of the X_{mix(pj)},
target_mix_x a data.frame summarizing the target distributions of the X_{mix(pj)}
If Poisson variables are produced:
pois_sum_x a data.frame summarizing the simulated distributions of the X_{pois(pj)}
If Negative Binomial variables are produced:
nb_sum_x a data.frame summarizing the simulated distributions of the X_{nb(pj)}
If random effects are produced:
cont_sum_u a data.frame summarizing the simulated distributions of the U_{cont(pj)} and U_{comp(pj)},
target_sum_u a data.frame summarizing the target distributions of the U_{cont(pj)} and U_{comp(pj)},
sum_uall a data.frame summarizing the simulated distributions of U_all,
mix_sum_u a data.frame summarizing the simulated distributions of the U_{mix(pj)},
target_mix_u a data.frame summarizing the target distributions of the U_{mix(pj)},
rho.u list of length M, each component a list of length M;
rho.u[[p]][[q]] = cor(cbind(U[[p]], U[[q]])) if p != q or rho.u[[p]][[q]] = cor(U[[p]])) if p = q
rho.uall list of length M, each component a list of length M;
rho.uall[[p]][[q]] = cor(cbind(U_all[[p]], U_all[[q]])) if p != q or rho.uall[[p]][[q]] = cor(U_all[[p]])) if p = q
maxerr_u list of length M containing a vector of length M with the maximum correlation errors for U between outcomes
maxerr_u[[p]]][q] = abs(max(corr.u[[p]][[q]] - rho.u[[p]][[q]]))
References
See references for SimRepeat.
See Also
Examples
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 | M <- 3
B <- calc_theory("Beta", c(4, 1.5))
skews <- lapply(seq_len(M), function(x) B[3])
skurts <- lapply(seq_len(M), function(x) B[4])
fifths <- lapply(seq_len(M), function(x) B[5])
sixths <- lapply(seq_len(M), function(x) B[6])
Six <- lapply(seq_len(M), function(x) list(0.03))
corr.e <- matrix(c(1, 0.4, 0.4^2, 0.4, 1, 0.4, 0.4^2, 0.4, 1), M, M,
byrow = TRUE)
means <- lapply(seq_len(M), function(x) B[1])
vars <- lapply(seq_len(M), function(x) B[2]^2)
marginal <- list(0.3, 0.4, 0.5)
support <- lapply(seq_len(M), function(x) list(0:1))
corr.x <- list(list(matrix(1, 1, 1), matrix(0.4, 1, 1), matrix(0.4, 1, 1)),
list(matrix(0.4, 1, 1), matrix(1, 1, 1), matrix(0.4, 1, 1)),
list(matrix(0.4, 1, 1), matrix(0.4, 1, 1), matrix(1, 1, 1)))
betas <- list(0.5)
betas.t <- 1
betas.tint <- list(0.25)
Sys1 <- corrsys(10000, M, Time = 1:M, "Polynomial", "non_mix", means, vars,
skews, skurts, fifths, sixths, Six, marginal = marginal, support = support,
corr.x = corr.x, corr.e = corr.e, betas = betas, betas.t = betas.t,
betas.tint = betas.tint, quiet = TRUE)
Sum1 <- summary_sys(Sys1$Y, Sys1$E, E_mix = NULL, Sys1$X, Sys1$X_all, M,
"Polynomial", means, vars, skews, skurts, fifths, sixths,
marginal = marginal, support = support, corr.x = corr.x, corr.e = corr.e)
## Not run:
seed <- 276
n <- 10000
M <- 3
Time <- 1:M
# Error terms have a beta(4, 1.5) distribution with an AR(1, p = 0.4)
# correlation structure
B <- calc_theory("Beta", c(4, 1.5))
skews <- lapply(seq_len(M), function(x) B[3])
skurts <- lapply(seq_len(M), function(x) B[4])
fifths <- lapply(seq_len(M), function(x) B[5])
sixths <- lapply(seq_len(M), function(x) B[6])
Six <- lapply(seq_len(M), function(x) list(0.03))
error_type <- "non_mix"
corr.e <- matrix(c(1, 0.4, 0.4^2, 0.4, 1, 0.4, 0.4^2, 0.4, 1), M, M,
byrow = TRUE)
# 1 continuous mixture of Normal(-2, 1) and Normal(2, 1) for each Y
mix_pis <- lapply(seq_len(M), function(x) list(c(0.4, 0.6)))
mix_mus <- lapply(seq_len(M), function(x) list(c(-2, 2)))
mix_sigmas <- lapply(seq_len(M), function(x) list(c(1, 1)))
mix_skews <- lapply(seq_len(M), function(x) list(c(0, 0)))
mix_skurts <- lapply(seq_len(M), function(x) list(c(0, 0)))
mix_fifths <- lapply(seq_len(M), function(x) list(c(0, 0)))
mix_sixths <- lapply(seq_len(M), function(x) list(c(0, 0)))
mix_Six <- list()
Nstcum <- calc_mixmoments(mix_pis[[1]][[1]], mix_mus[[1]][[1]],
mix_sigmas[[1]][[1]], mix_skews[[1]][[1]], mix_skurts[[1]][[1]],
mix_fifths[[1]][[1]], mix_sixths[[1]][[1]])
means <- lapply(seq_len(M), function(x) c(Nstcum[1], B[1]))
vars <- lapply(seq_len(M), function(x) c(Nstcum[2]^2, B[2]^2))
# 1 binary variable for each Y
marginal <- lapply(seq_len(M), function(x) list(0.4))
support <- list(NULL, list(c(0, 1)), NULL)
# 1 Poisson variable for each Y
lam <- list(1, 5, 10)
# Y2 and Y3 are zero-inflated Poisson variables
p_zip <- list(NULL, 0.05, 0.1)
# 1 NB variable for each Y
size <- list(10, 15, 20)
prob <- list(0.3, 0.4, 0.5)
# either prob or mu is required (not both)
mu <- mapply(function(x, y) x * (1 - y)/y, size, prob, SIMPLIFY = FALSE)
# Y2 and Y3 are zero-inflated NB variables
p_zinb <- list(NULL, 0.05, 0.1)
# The 2nd (the normal mixture) variable is the same across Y
same.var <- 2
# Create the correlation matrix in terms of the components of the normal
# mixture
K <- 5
corr.x <- list()
corr.x[[1]] <- list(matrix(0.1, K, K), matrix(0.2, K, K), matrix(0.3, K, K))
diag(corr.x[[1]][[1]]) <- 1
# set correlation between components to 0
corr.x[[1]][[1]][2:3, 2:3] <- diag(2)
# set correlations with the same variable equal across outcomes
corr.x[[1]][[2]][, same.var] <- corr.x[[1]][[3]][, same.var] <-
corr.x[[1]][[1]][, same.var]
corr.x[[2]] <- list(t(corr.x[[1]][[2]]), matrix(0.35, K, K),
matrix(0.4, K, K))
diag(corr.x[[2]][[2]]) <- 1
corr.x[[2]][[2]][2:3, 2:3] <- diag(2)
corr.x[[2]][[2]][, same.var] <- corr.x[[2]][[3]][, same.var] <-
t(corr.x[[1]][[2]][same.var, ])
corr.x[[2]][[3]][same.var, ] <- corr.x[[1]][[3]][same.var, ]
corr.x[[2]][[2]][same.var, ] <- t(corr.x[[2]][[2]][, same.var])
corr.x[[3]] <- list(t(corr.x[[1]][[3]]), t(corr.x[[2]][[3]]),
matrix(0.5, K, K))
diag(corr.x[[3]][[3]]) <- 1
corr.x[[3]][[3]][2:3, 2:3] <- diag(2)
corr.x[[3]][[3]][, same.var] <- t(corr.x[[1]][[3]][same.var, ])
corr.x[[3]][[3]][same.var, ] <- t(corr.x[[3]][[3]][, same.var])
# The 2nd and 3rd variables of each Y are subject-level variables
subj.var <- matrix(c(1, 2, 1, 3, 2, 2, 2, 3, 3, 2, 3, 3), 6, 2, byrow = TRUE)
int.var <- tint.var <- NULL
betas.0 <- 0
betas <- list(seq(0.5, 0.5 + (K - 2) * 0.25, 0.25))
betas.subj <- list(seq(0.5, 0.5 + (K - 2) * 0.1, 0.1))
betas.int <- list()
betas.t <- 1
betas.tint <- list(c(0.25, 0.5))
method <- "Polynomial"
# Check parameter inputs
checkpar(M, method, error_type, means, vars, skews, skurts, fifths, sixths,
Six, mix_pis, mix_mus, mix_sigmas, mix_skews, mix_skurts, mix_fifths,
mix_sixths, mix_Six, marginal, support, lam, p_zip, pois_eps = list(),
size, prob, mu, p_zinb, nb_eps = list(), corr.x, corr.yx = list(),
corr.e, same.var, subj.var, int.var, tint.var, betas.0, betas,
betas.subj, betas.int, betas.t, betas.tint)
# Simulated system using correlation method 1
N <- corrsys(n, M, Time, method, error_type, means, vars, skews, skurts,
fifths, sixths, Six, mix_pis, mix_mus, mix_sigmas, mix_skews, mix_skurts,
mix_fifths, mix_sixths, mix_Six, marginal, support, lam, p_zip, size,
prob, mu, p_zinb, corr.x, corr.e, same.var, subj.var, int.var, tint.var,
betas.0, betas, betas.subj, betas.int, betas.t, betas.tint, seed = seed,
use.nearPD = FALSE)
# Summarize the results
S <- summary_sys(N$Y, N$E, E_mix = NULL, N$X, N$X_all, M, method, means,
vars, skews, skurts, fifths, sixths, mix_pis, mix_mus, mix_sigmas,
mix_skews, mix_skurts, mix_fifths, mix_sixths, marginal, support, lam,
p_zip, size, prob, mu, p_zinb, corr.x, corr.e)
S$sum_xall
S$maxerr
## End(Not run)
|
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.