Skip to content

Commit ada8899

Browse files
ran stylr through it
1 parent 117679d commit ada8899

File tree

12 files changed

+435
-420
lines changed

12 files changed

+435
-420
lines changed

DESCRIPTION

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -8,7 +8,7 @@ Description: A few functions aim to provide a statistic tool for three purposes.
88
License: MIT + file LICENSE
99
Encoding: UTF-8
1010
Roxygen: list(markdown = TRUE)
11-
RoxygenNote: 7.1.2
11+
RoxygenNote: 7.3.2
1212
Imports:
1313
OpenMx (>= 2.19.6),
1414
stats (>= 3.5.0)

R/Internal.R

Lines changed: 8 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -4,9 +4,12 @@
44
#' @param sigma Covariance matrix
55
#' @return Generates multivariate normal data from a covariance matrix (\code{sigma}) of length \code{n}
66
#'
7-
rmvn <- function(n,sigma) {
8-
Sh <- with(svd(sigma),
9-
v%*%diag(sqrt(d))%*%t(u))
10-
matrix(stats::rnorm(ncol(sigma)*n),
11-
ncol = ncol(sigma))%*%Sh
7+
rmvn <- function(n, sigma) {
8+
Sh <- with(
9+
svd(sigma),
10+
v %*% diag(sqrt(d)) %*% t(u)
11+
)
12+
matrix(stats::rnorm(ncol(sigma) * n),
13+
ncol = ncol(sigma)
14+
) %*% Sh
1215
}

R/Power_LS.R

Lines changed: 25 additions & 30 deletions
Original file line numberDiff line numberDiff line change
@@ -17,34 +17,29 @@
1717
#' @export
1818

1919

20-
Power_LS <- function(N1,N2,power,p_N1=NULL,h2,c2,R1=1,R2=.5,alpha = .05){
21-
Za <- qnorm(1-alpha)
22-
23-
if(missing(power)){
24-
Zb <- sqrt(h2^2*(abs(R1-R2)^2) / ((1-(R1*h2+c2)^2)^2/N1 + (1-(R2*h2+c2)^2)^2/N2)) - Za
25-
power_result <- pnorm(Zb,0)
26-
return(power_result)
27-
}
28-
if(!missing(N1) & missing(N2)){
29-
Zb <- qnorm(power)
30-
N2_result <- (1-(R2*h2+c2)^2)^2 / (h2^2*(abs(R1-R2)^2)/((Za+Zb)^2) - ((1-(R1*h2+c2)^2)^2/N1))
31-
return(round(N2_result))
32-
}
33-
if(missing(N1) & !missing(N2)){
34-
Zb <- qnorm(power)
35-
N1_result <- (1-(R1*h2+c2)^2)^2 / (h2^2*(abs(R1-R2)^2)/((Za+Zb)^2) - ((1-(R2*h2+c2)^2)^2/N2))
36-
return(round(N1_result))
37-
}
38-
if(missing(N1) & missing(N2) & !is.null(p_N1)){
39-
Zb <- qnorm(power)
40-
N_total <- (1-(R1*h2+c2)^2)^2/(p_N1*h2^2*abs(R1-R2)^2/((Za+Zb)^2)) + (1-(R2*h2+c2)^2)^2/((1-p_N1)*h2^2*(abs(R1-R2)^2)/((Za+Zb)^2))
41-
N1_result <- N_total * p_N1
42-
N2_result <- N_total * (1-p_N1)
43-
return(c(round(N1_result),round(N2_result)))
44-
45-
}
46-
47-
20+
Power_LS <- function(N1, N2, power, p_N1 = NULL, h2, c2, R1 = 1, R2 = .5, alpha = .05) {
21+
Za <- qnorm(1 - alpha)
22+
23+
if (missing(power)) {
24+
Zb <- sqrt(h2^2 * (abs(R1 - R2)^2) / ((1 - (R1 * h2 + c2)^2)^2 / N1 + (1 - (R2 * h2 + c2)^2)^2 / N2)) - Za
25+
power_result <- pnorm(Zb, 0)
26+
return(power_result)
27+
}
28+
if (!missing(N1) & missing(N2)) {
29+
Zb <- qnorm(power)
30+
N2_result <- (1 - (R2 * h2 + c2)^2)^2 / (h2^2 * (abs(R1 - R2)^2) / ((Za + Zb)^2) - ((1 - (R1 * h2 + c2)^2)^2 / N1))
31+
return(round(N2_result))
32+
}
33+
if (missing(N1) & !missing(N2)) {
34+
Zb <- qnorm(power)
35+
N1_result <- (1 - (R1 * h2 + c2)^2)^2 / (h2^2 * (abs(R1 - R2)^2) / ((Za + Zb)^2) - ((1 - (R2 * h2 + c2)^2)^2 / N2))
36+
return(round(N1_result))
37+
}
38+
if (missing(N1) & missing(N2) & !is.null(p_N1)) {
39+
Zb <- qnorm(power)
40+
N_total <- (1 - (R1 * h2 + c2)^2)^2 / (p_N1 * h2^2 * abs(R1 - R2)^2 / ((Za + Zb)^2)) + (1 - (R2 * h2 + c2)^2)^2 / ((1 - p_N1) * h2^2 * (abs(R1 - R2)^2) / ((Za + Zb)^2))
41+
N1_result <- N_total * p_N1
42+
N2_result <- N_total * (1 - p_N1)
43+
return(c(round(N1_result), round(N2_result)))
44+
}
4845
}
49-
50-

R/Sim_Fit.R

Lines changed: 47 additions & 44 deletions
Original file line numberDiff line numberDiff line change
@@ -17,48 +17,51 @@
1717
#' \item{Data}{A \code{data.frame} consists of the simulated raw data}
1818
#' @export
1919

20-
Sim_Fit <- function(
21-
GroupNames = c("KinPair1","KinPair2"),
22-
GroupSizes = c(100,100),
23-
nIter = 100,
24-
SSeed = 62,
25-
GroupRel=c(1,.5),
26-
GroupR_c=c(1,1),
27-
mu = c(0,0),
28-
ace1=c(1,1,1),
29-
ace2=c(1,1,1),
30-
ifComb=FALSE,
31-
lbound=FALSE,
32-
saveRaw=FALSE
33-
){
34-
l.results <- list()
35-
for(i in 1: nIter){
36-
set.seed(SSeed - 1 + i)
37-
df_temp <- kinsim_double(
38-
GroupNames = GroupNames,
39-
GroupSizes = GroupSizes,
40-
GroupRel=GroupRel,
41-
GroupR_c=GroupR_c,
42-
mu = mu,
43-
ace1=ace1,
44-
ace2=ace2,
45-
ifComb=ifComb)
46-
if(!saveRaw){
47-
l.results[[i]] <- list(
48-
Results = fit_uniACE(
49-
data_1 = df_temp[which(df_temp$GroupName == GroupNames[1]), c("y1","y2")],
50-
data_2 = df_temp[which(df_temp$GroupName == GroupNames[2]), c("y1","y2")],
51-
GroupRel = GroupRel, GroupR_c = GroupR_c, lbound = lbound),
52-
data = NA)
53-
}else{
54-
l.results[[i]] <- list(
55-
Results = fit_uniACE(
56-
data_1 = df_temp[which(df_temp$GroupName == GroupNames[1]), c("y1","y2")],
57-
data_2 = df_temp[which(df_temp$GroupName == GroupNames[2]), c("y1","y2")],
58-
GroupRel = GroupRel, GroupR_c = GroupR_c, lbound = lbound),
59-
data = df_temp)
60-
}
61-
names(l.results)[[i]] <- paste("Iteration",i, sep = "")
62-
}
63-
return(l.results)
20+
Sim_Fit <- function(GroupNames = c("KinPair1", "KinPair2"),
21+
GroupSizes = c(100, 100),
22+
nIter = 100,
23+
SSeed = 62,
24+
GroupRel = c(1, .5),
25+
GroupR_c = c(1, 1),
26+
mu = c(0, 0),
27+
ace1 = c(1, 1, 1),
28+
ace2 = c(1, 1, 1),
29+
ifComb = FALSE,
30+
lbound = FALSE,
31+
saveRaw = FALSE) {
32+
l.results <- list()
33+
for (i in 1:nIter) {
34+
set.seed(SSeed - 1 + i)
35+
df_temp <- kinsim_double(
36+
GroupNames = GroupNames,
37+
GroupSizes = GroupSizes,
38+
GroupRel = GroupRel,
39+
GroupR_c = GroupR_c,
40+
mu = mu,
41+
ace1 = ace1,
42+
ace2 = ace2,
43+
ifComb = ifComb
44+
)
45+
if (!saveRaw) {
46+
l.results[[i]] <- list(
47+
Results = fit_uniACE(
48+
data_1 = df_temp[which(df_temp$GroupName == GroupNames[1]), c("y1", "y2")],
49+
data_2 = df_temp[which(df_temp$GroupName == GroupNames[2]), c("y1", "y2")],
50+
GroupRel = GroupRel, GroupR_c = GroupR_c, lbound = lbound
51+
),
52+
data = NA
53+
)
54+
} else {
55+
l.results[[i]] <- list(
56+
Results = fit_uniACE(
57+
data_1 = df_temp[which(df_temp$GroupName == GroupNames[1]), c("y1", "y2")],
58+
data_2 = df_temp[which(df_temp$GroupName == GroupNames[2]), c("y1", "y2")],
59+
GroupRel = GroupRel, GroupR_c = GroupR_c, lbound = lbound
60+
),
61+
data = df_temp
62+
)
63+
}
64+
names(l.results)[[i]] <- paste("Iteration", i, sep = "")
65+
}
66+
return(l.results)
6467
}

0 commit comments

Comments
 (0)