From d62f7a8b63f3033d870f5b24fee5abb4f5df3b36 Mon Sep 17 00:00:00 2001
From: Kaylea Champion <khascall@gmail.com>
Date: Tue, 14 Nov 2023 16:12:51 -0800
Subject: [PATCH] stubs in a place to insert the effect size

---
 R/powerAnalysis.R | 23 +++++++++++++++--------
 1 file changed, 15 insertions(+), 8 deletions(-)

diff --git a/R/powerAnalysis.R b/R/powerAnalysis.R
index 724205c..8eb0c44 100644
--- a/R/powerAnalysis.R
+++ b/R/powerAnalysis.R
@@ -17,7 +17,8 @@ l2p <- function(b) {
 #Matt: 
 makeDataNew <- function(n) {
   tDF <- data.frame(
-    up.fac.mean=rnorm(n=n, mean=-0.1296376, sd=1.479847), # up.fac.mean
+    ## don't sim the outcome
+    #up.fac.mean=rnorm(n=n, mean=-0.1296376, sd=1.479847), # up.fac.mean
     mmt=rlnorm(n=n, mean=1.685715, sd = 0.2532059), # mmt
     ## this generates a 50-50 split of milestones --v
     #milestones=rbinom(n=n, size=1, prob=c(0.247, 0.753)), #milestones
@@ -32,7 +33,8 @@ makeDataNew <- function(n) {
 
 makeDataNew2 <- function(n) {
   tDF <- data.frame(
-    up.fac.mean=rnorm(n=n, mean=-0.1296376, sd=1.479847), # up.fac.mean
+    ## don't sim the outcome
+    #up.fac.mean=rnorm(n=n, mean=-0.1296376, sd=1.479847), # up.fac.mean
     formal.score=rlnorm(n=n, mean=6.220282, sd = 2.544058) # formal.score
   )
   tDF[is.na(tDF) | tDF=="Inf"] = NA
@@ -51,18 +53,22 @@ powerCheck <- function(n, nSims) { #run a power calculation on the dataset given
   signifM <- rep(NA, nSims)
   for (s in 1:nSims) {           # repeatedly we will....
     simData <- makeDataNew(n)       # make some data
+    ## outcome goes here --v
+    # e.g. simData$up.fac.mean <- (usefuleffsizeA * mmt) + (usefuleffsizeB * milestones) + rnorm(n=1, mean=0, sd=1) ##plus some noise
     #have updated for kkex through here, now need to look at the underproduction work
     #m1.sim <- lm(up.fac.mean ~ ((mmt)/ (milestones/age)), data=simData)
-    m1.sim <- lm(up.fac.mean ~ mmt + milestones + age, data=simData)
+    ## could leave age out for now?
+    #m1.sim <- lm(up.fac.mean ~ mmt + milestones + age, data=simData)
+    m1.sim <- lm(up.fac.mean ~ mmt + milestones, data=simData)
     p0 <- coef(summary(m1.sim))[1,4] #intercept
     p1 <- coef(summary(m1.sim))[2,4] #mmt
     p2 <- coef(summary(m1.sim))[3,4] #milestones
-    p3 <- coef(summary(m1.sim))[4,4] #age
+    #p3 <- coef(summary(m1.sim))[4,4] #age
     signif0[s] <- p0 <=.05
     signif1[s] <- p1 <=.05
     signif2[s] <- p2 <=.05
-    signif3[s] <- p3 <=.05
-    signifM[s] <- p0 <=.05 & p1 <=.05 & p2 <=.05 & p3 <=.05
+    #signif3[s] <- p3 <=.05
+    signifM[s] <- p0 <=.05 & p1 <=.05 & p2 <=.05 #& p3 <=.05
   }  
   power <- c(mean(signif0), mean(signif1), mean(signif2), mean(signif3), mean(signifM))
   return(power)
@@ -77,14 +83,15 @@ powerCheck2 <- function(n, nSims) { #run a power calculation on the dataset give
     simData <- makeDataNew2(n)       # make some data
     #have updated for kkex through here, now need to look at the underproduction work
     #m1.sim <- lm(up.fac.mean ~ ((mmt)/ (milestones/age)), data=simData)
+    ## outcome goes here --v
+    ## e.g. simData$up.fac.mean <- (usefuleffsizeC * formal.score) + rnorm(1, mean=0, sd=1) ##plus some noise
     m1.sim <- lm(up.fac.mean ~ formal.score, data=simData)
     p0 <- coef(summary(m1.sim))[1,4] 
     p1 <- coef(summary(m1.sim))[2,4]
-    signif0[s] <- p0 <=.05
+    signif0[s] <- p0 <=.05 
     signif1[s] <- p1 <=.05
     signifM[s] <- p0 <=.05 & p1 <=.05 
   }  
   power <- c(mean(signif0), mean(signif1), mean(signifM))
   return(power)
 }
-