On this page,we take an example to evaluate the numerical performance of our proposed G-CRT (Algorithm 4) with the statistic functions discussed in Section 4.3 for testing conditional independence. We employ a linear regression model for \(\mathcal{L}(Y\mid X)\) with a low-dimension setting (\(p=20\) and \(n=50\)) and fix the procedure parameters \(M=100\) and \(L=3\) for conducting Algorithm 2.
library(glmnet)
library(rags2ridges)
library(glasso)
library(Matrix)
library(MASS)
library(mvtnorm)
library(reshape)
library(ggplot2)
source("CRT-Simulate.R")
source("CRT-Stat.R")
source("DBL_GLMs_functions.R")
source("OneExp.R")
source("Run-Testing.R")
source("dCRT_function.R")
source("GoF-Settings.R")
source("CSS-GGM.R")
source("GGMTesting.R")
source("CIT_function_for_webpage.R")In this setting, we examine our G-CRT procedure with the statistics LM-SST and LM-SSR by comparison to the classical F-test and the Bonferroni-adjusted dCRT using Gaussian Lasso models.
# Set parameters
popul_param_llw <- list(
setting = "l_l_w",
N = c(50),
P = c(20),
setT=1:8,
modelParam = list(
K = 6,
s = 0.2,
grid_signal= seq(0,2.5,by=0.25)
)
)Repeat experiments 10 times
path_prefix=paste0('Experiments-llw','/')
if(!dir.exists(path_prefix)){
dir.create(path_prefix)
}
# Repeat experiment and save results
for (epo in 1:10){
output = CRT_webpage(popul_param_llw,epo = epo)
save(output,popul_param_llw,file=paste0(path_prefix,epo,'.RData'))
}
# Load results
all_output=list()
for(epo in 1:10){
load(file=paste0(path_prefix,epo,'.RData'))
all_output[[epo]]=output
}
combined_pvalues=list()
for(ind.signal in 1:length(popul_param_llw$modelParam$grid_signal)){
combined_pvalues[[ind.signal]]=t(sapply(all_output, function(v)v[[ind.signal]]$pvalues))
}
combined_pvalues=lapply(combined_pvalues,function(x)x[,-c(5)])
combined_pvalues=lapply(combined_pvalues,function(x){
colnames(x)=c("LM-SST","LM-SSR","F-test","dCRT")
return(x)})Present the power and standard error of various testing methods in the experiment, in which power is presented against varying values of \(\theta\)
Fun_SE=function(x){sd(x)/sqrt(length(x))}
pvalues_array=simplify2array(combined_pvalues)
sig.lv=0.05
power_table=data.frame(apply(pvalues_array<sig.lv, 3:2,mean));
power_SE_table=data.frame(apply(pvalues_array<sig.lv, 3:2,Fun_SE))
num_methods=ncol(power_table)
power_table$theta=popul_param_llw$modelParam$grid_signal
power_SE_table$theta=popul_param_llw$modelParam$grid_signal
df_melt=melt(power_table,id.vars="theta", variable_name = "Test")
df_melt$Power=df_melt$value;df_melt$value=NULL
SE_melt <- melt(power_SE_table, id.vars = "theta", variable_name = "Test")
SE_melt$SE=SE_melt$value;SE_melt$value=NULL
# Merge the data and SE data frames
merged_df <- merge(df_melt, SE_melt, by = c("theta", "Test"))
print(merged_df)## theta Test Power SE
## 1 0.00 dCRT 0.1 0.1000000
## 2 0.00 F.test 0.0 0.0000000
## 3 0.00 LM.SSR 0.0 0.0000000
## 4 0.00 LM.SST 0.0 0.0000000
## 5 0.25 dCRT 0.1 0.1000000
## 6 0.25 F.test 0.0 0.0000000
## 7 0.25 LM.SSR 0.0 0.0000000
## 8 0.25 LM.SST 0.0 0.0000000
## 9 0.50 dCRT 0.1 0.1000000
## 10 0.50 F.test 0.2 0.1333333
## 11 0.50 LM.SSR 0.2 0.1333333
## 12 0.50 LM.SST 0.1 0.1000000
## 13 0.75 dCRT 0.0 0.0000000
## 14 0.75 F.test 0.3 0.1527525
## 15 0.75 LM.SSR 0.3 0.1527525
## 16 0.75 LM.SST 0.3 0.1527525
## 17 1.00 dCRT 0.2 0.1333333
## 18 1.00 F.test 0.8 0.1333333
## 19 1.00 LM.SSR 0.8 0.1333333
## 20 1.00 LM.SST 0.9 0.1000000
## 21 1.25 dCRT 0.2 0.1333333
## 22 1.25 F.test 1.0 0.0000000
## 23 1.25 LM.SSR 1.0 0.0000000
## 24 1.25 LM.SST 1.0 0.0000000
## 25 1.50 dCRT 0.4 0.1632993
## 26 1.50 F.test 1.0 0.0000000
## 27 1.50 LM.SSR 1.0 0.0000000
## 28 1.50 LM.SST 1.0 0.0000000
## 29 1.75 dCRT 0.4 0.1632993
## 30 1.75 F.test 1.0 0.0000000
## 31 1.75 LM.SSR 1.0 0.0000000
## 32 1.75 LM.SST 1.0 0.0000000
## 33 2.00 dCRT 0.7 0.1527525
## 34 2.00 F.test 1.0 0.0000000
## 35 2.00 LM.SSR 1.0 0.0000000
## 36 2.00 LM.SST 1.0 0.0000000
## 37 2.25 dCRT 0.6 0.1632993
## 38 2.25 F.test 1.0 0.0000000
## 39 2.25 LM.SSR 1.0 0.0000000
## 40 2.25 LM.SST 1.0 0.0000000
## 41 2.50 dCRT 0.9 0.1000000
## 42 2.50 F.test 1.0 0.0000000
## 43 2.50 LM.SSR 1.0 0.0000000
## 44 2.50 LM.SST 1.0 0.0000000