sim.mpVSp <- function(samples=500, nobs=100, time = 4, nvar = 2, AB = 5, nn = 10, type = "rook", beta = 1 ){ source("predictors.R") source("kpt.panel.R") library(splm) library(sphet) rhos<- c(-0.4, 0.4) lambdas<- c(-0.6, -0.4, -0.2, 0.2, 0.4, 0.6) nobs <- nobs time <- time nvar <- nvar NT <- nobs*time coord <- cbind(rep(seq(1, sqrt(nobs)), each = sqrt(nobs)), rep(seq(1,sqrt(nobs)),sqrt(nobs))) Wa <- listw2dgCMatrix(nb2listw(circular(sqrt(nobs),sqrt(nobs), AB))) Wb <- listw2dgCMatrix(nb2listw(knn2nb(knearneigh(coord, k = nn)))) Wc <- listw2dgCMatrix(nb2listw(cell2nb(sqrt(nobs),sqrt(nobs), type = type))) Itime <- Diagonal(time) IN <- Diagonal(nobs) INT <- Diagonal(NT) IWa <- kronecker(Itime, Wa) IWb <- kronecker(Itime, Wb) IWc <- kronecker(Itime, Wc) beta <- rep(beta, nvar) mu11 <- runif(time,0,1) mu12 <- runif(time,0,1) x11 <- vector("numeric", length = 0) x12 <- vector("numeric", length = 0) for (i in 1:time) x11 <- c(x11, rnorm(nobs, mu11[i], 1)) for (i in 1:time) x12 <- c(x12, rnorm(nobs, mu12[i], 1)) xtot <- cbind(x11, x12) xtotb <- xtot %*% beta Rtilde <- var(xtotb)/(var(xtotb)+1) print(Rtilde) mu21 <- runif(time,0,1) x21 <- vector("numeric", length = 0) for (i in 1:time) x21 <- c(x21, rnorm(nobs, mu21[i], 1)) x22 <- .5 * x12 + rnorm(NT, 0, 1) xtot1 <- cbind(x21, x22) # cov(x12,x22)/(sqrt(var(x12))*sqrt(var(x22))) xtot1b <- xtot1 %*% beta indic <- rep(1:nobs,time) tindic <- rep(1:time, each = nobs) index <- cbind(indic, tindic) Xwithin <- splm:::panel.transformations(xtot,indic, type= "within") X1within <- splm:::panel.transformations(xtot1,indic, type= "within") #1 WAXwithin <- IWa %*% Xwithin WWAXwithin <- IWa %*% WAXwithin WAX1within <- IWa %*% X1within WWAX1within <- IWa %*% WAX1within #2 WBXwithin <- IWb %*% Xwithin WWBXwithin <- IWb %*% WBXwithin WBX1within <- IWb %*% X1within WWBX1within <- IWb %*% WBX1within #3 WCXwithin <- IWc %*% Xwithin WWCXwithin <- IWc %*% WCXwithin WCX1within <- IWc %*% X1within WWCX1within <- IWc %*% WCX1within Hexp1 <- cbind(as.matrix(Xwithin), as.matrix(WAXwithin), as.matrix(WBXwithin), as.matrix(WWAXwithin), as.matrix(WWBXwithin)) Hexp2 <- cbind(as.matrix(Xwithin), as.matrix(WAXwithin), as.matrix(WCXwithin), as.matrix(WWAXwithin), as.matrix(WWCXwithin)) Hexp3 <- cbind(as.matrix(Xwithin), as.matrix(X1within), as.matrix(WAXwithin), as.matrix(WAX1within), as.matrix(WWAXwithin),as.matrix(WWAX1within)) Hexp4 <- cbind(as.matrix(Xwithin), as.matrix(X1within), as.matrix(WAXwithin), as.matrix(WCX1within), as.matrix(WWAXwithin), as.matrix(WWCX1within)) Hexp5 <- cbind(as.matrix(Xwithin), as.matrix(WAXwithin), as.matrix(WBXwithin), as.matrix(WCXwithin), as.matrix(WWAXwithin), as.matrix(WWBXwithin), as.matrix(WWCXwithin)) Hexp6 <- cbind(as.matrix(Xwithin), as.matrix(X1within), as.matrix(WAXwithin), as.matrix(WBX1within), as.matrix(WCX1within), as.matrix(WWAXwithin), as.matrix(WWBX1within), as.matrix(WWCX1within)) # ################################################# ### First Set of Experiments ################################################# FIN1<-vector("numeric",length=0) FIN2<-vector("numeric",length=0) FIN3<-vector("numeric",length=0) FIN4<-vector("numeric",length=0) FIN5<-vector("numeric",length=0) FIN6<-vector("numeric",length=0) FIN7<-vector("numeric",length=0) FIN8<-vector("numeric",length=0) for(rho in rhos){ IrWA<-INT - rho*IWa IrWB<-INT - rho*IWb IrWAinv<-solve(IrWA) IrWBinv<-solve(IrWB) TMP1<-vector("numeric",length=0) TMP2<-vector("numeric",length=0) TMP3<-vector("numeric",length=0) TMP4<-vector("numeric",length=0) TMP5<-vector("numeric",length=0) TMP6<-vector("numeric",length=0) TMP7<-vector("numeric",length=0) TMP8<-vector("numeric",length=0) for (lambda in lambdas){ IlWA<-INT - lambda*IWa IlWB<-INT - lambda*IWb IlWAinv<-solve(IlWA) IlWBinv<-solve(IlWB) tmp1<-vector("numeric",length=samples) tmp2<-vector("numeric",length=samples) tmp3<-vector("numeric",length=samples) tmp4<-vector("numeric",length=samples) tmp5<-vector("numeric",length=samples) tmp6<-vector("numeric",length=samples) tmp7<-vector("numeric",length=samples) tmp8<-vector("numeric",length=samples) i=1 while(i<=samples){ vit0 <- rnorm(NT,mean = 0, sd = 1) vit1 <- rnorm(NT,mean = 0, sd = 1) vit0nn <- (exp(rnorm(NT,0,1)) - exp(0.5)) /sqrt(exp(2)-exp(1)) vit1nn <- (exp(rnorm(NT,0,1)) - exp(0.5)) /sqrt(exp(2)-exp(1)) uit0 <- IrWAinv %*% vit0 uit1 <- vit1 uit0nn <- IrWAinv %*% vit0nn uit1nn <- vit1nn y0<- as.matrix(IlWAinv %*% (xtotb + uit0)) y1<- as.matrix(IlWBinv %*% (xtotb + uit1)) y0nn <- as.matrix(IlWAinv %*% (xtotb + uit0nn)) y1nn <- as.matrix(IlWBinv %*% (xtotb + uit1nn)) colnames(y0) <- "y0" colnames(y1) <- "y1" colnames(y0nn) <- "y0nn" colnames(y1nn) <- "y1nn" data1 <- data.frame(cbind(index, y0, y1, y0nn, y1nn, xtot)) H0_size <- spgm(y0 ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T) H1_size <- spgm(y0 ~ x11 + x12, data = data1, index = NULL, listw = Wb, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H0_power <- spgm(y1 ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T) H1_power <- spgm(y1 ~ x11 + x12, data = data1, index = NULL, listw = Wb, model = "within", lag = T, spatial.error =F, method = "w2sls") H0_sizenn <- spgm(y0nn ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T) H1_sizenn <- spgm(y0nn ~ x11 + x12, data = data1, index = NULL, listw = Wb, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H0_powernn <- spgm(y1nn ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T) H1_powernn <- spgm(y1nn ~ x11 + x12, data = data1, index = NULL, listw = Wb, model = "within", lag = T, spatial.error =F, method = "w2sls") pred_size<-predictors(H0_size, H1_size, y = y0, xreg = xtot, w0 = IWa, w1 = IWb, nobs = nobs, time = time, IN = IN, IT = Itime) pred_power<-predictors(H0_power, H1_power, y = y1, xreg = xtot, w0 = IWa, w1 = IWb, nobs = nobs, time = time, IN = IN, IT = Itime) pred_sizenn<-predictors(H0_sizenn, H1_sizenn, y = y0nn, xreg = xtot, w0 = IWa, w1 = IWb, nobs = nobs, time = time, IN = IN, IT = Itime) pred_powernn<-predictors(H0_powernn, H1_powernn, y = y1nn, xreg = xtot, w0 = IWa, w1 = IWb, nobs = nobs, time = time, IN = IN, IT = Itime) test_size_p1<-try(kptest(y0,H0_size, pred = pred_size[,1], x0 = xtot, w0 = IWa, instruments = Hexp1, nobs = nobs, time = time)) test_size_p2<-try(kptest(y0, H0_size, pred = pred_size[,2], x0 = xtot, w0 = IWa, instruments = Hexp1, nobs = nobs, time = time)) test_power_p1<-try(kptest(y1,H0_power, pred=pred_power[,1], x0 = xtot, w0 = IWa, instruments = Hexp1, nobs = nobs, time = time)) test_power_p2<-try(kptest(y1,H0_power, pred=pred_power[,2], x0 = xtot, w0 = IWa, instruments = Hexp1, nobs = nobs, time = time)) test_size_p1nn<-try(kptest(y0nn,H0_sizenn, pred = pred_sizenn[,1], x0 = xtot, w0 = IWa, instruments = Hexp1, nobs = nobs, time = time)) test_size_p2nn<-try(kptest(y0nn, H0_sizenn, pred = pred_sizenn[,2], x0 = xtot, w0 = IWa, instruments = Hexp1, nobs = nobs, time = time)) test_power_p1nn<-try(kptest(y1nn,H0_powernn, pred=pred_powernn[,1], x0 = xtot, w0 = IWa, instruments = Hexp1, nobs = nobs, time = time)) test_power_p2nn<-try(kptest(y1nn,H0_powernn, pred=pred_powernn[,2], x0 = xtot, w0 = IWa, instruments = Hexp1, nobs = nobs, time = time)) tmp1[i]<-test_size_p1 tmp2[i]<-test_size_p2 tmp3[i]<-test_power_p1 tmp4[i]<-test_power_p2 tmp5[i]<-test_size_p1nn tmp6[i]<-test_size_p2nn tmp7[i]<-test_power_p1nn tmp8[i]<-test_power_p2nn i=i+1 if (i==samples) print(i) } TMP1<-c(TMP1, tmp1) TMP2<-c(TMP2, tmp2) TMP3<-c(TMP3, tmp3) TMP4<-c(TMP4, tmp4) TMP5<-c(TMP5, tmp5) TMP6<-c(TMP6, tmp6) TMP7<-c(TMP7, tmp7) TMP8<-c(TMP8, tmp8) } FIN1<-c(FIN1, TMP1) FIN2<-c(FIN2, TMP2) FIN3<-c(FIN3, TMP3) FIN4<-c(FIN4, TMP4) FIN5<-c(FIN5, TMP5) FIN6<-c(FIN6, TMP6) FIN7<-c(FIN7, TMP7) FIN8<-c(FIN8, TMP8) } STAT<- qchisq(0.05,1,lower.tail=FALSE) for (i in 1:length(FIN1)){ if(FIN1[i]>STAT) FIN1[i]=1 else FIN1[i]=0 if(FIN2[i]>STAT) FIN2[i]=1 else FIN2[i]=0 if(FIN3[i]>STAT) FIN3[i]=1 else FIN3[i]=0 if(FIN4[i]>STAT) FIN4[i]=1 else FIN4[i]=0 if(FIN5[i]>STAT) FIN5[i]=1 else FIN5[i]=0 if(FIN6[i]>STAT) FIN6[i]=1 else FIN6[i]=0 if(FIN7[i]>STAT) FIN7[i]=1 else FIN7[i]=0 if(FIN8[i]>STAT) FIN8[i]=1 else FIN8[i]=0 } org1<-matrix(FIN1,samples,length(lambdas)*length(rhos)) org2<-matrix(FIN2,samples,length(lambdas)*length(rhos)) org3<-matrix(FIN3,samples,length(lambdas)*length(rhos)) org4<-matrix(FIN4,samples,length(lambdas)*length(rhos)) org5<-matrix(FIN5,samples,length(lambdas)*length(rhos)) org6<-matrix(FIN6,samples,length(lambdas)*length(rhos)) org7<-matrix(FIN7,samples,length(lambdas)*length(rhos)) org8<-matrix(FIN8,samples,length(lambdas)*length(rhos)) FR1<-colSums(org1)/samples FR2<-colSums(org2)/samples FR3<-colSums(org3)/samples FR4<-colSums(org4)/samples FR5<-colSums(org5)/samples FR6<-colSums(org6)/samples FR7<-colSums(org7)/samples FR8<-colSums(org8)/samples table<-cbind(FR1,FR2, FR3, FR4, FR5,FR6, FR7, FR8) print(table) save(FIN1,FIN2,FIN3, FIN4, FIN5, FIN6, FIN7, FIN8, file = paste("exp1_",nobs,"_R2=",format(Rtilde, digits =2),sep="")) write.table(format(table, digits = 3, trim =TRUE),paste("tab1_",nobs,"_R2=",format(Rtilde, digits =2),".csv",sep=""), quote = FALSE, sep="&") ################################# ### Second Set of Experiments ###### ################################# FIN1<-vector("numeric",length=0) FIN2<-vector("numeric",length=0) FIN3<-vector("numeric",length=0) FIN4<-vector("numeric",length=0) FIN5<-vector("numeric",length=0) FIN6<-vector("numeric",length=0) FIN7<-vector("numeric",length=0) FIN8<-vector("numeric",length=0) for(rho in rhos){ IrWC<-INT - rho*IWc IrWA<-INT - rho*IWa IrWCinv<-solve(IrWC) IrWAinv<-solve(IrWA) TMP1<-vector("numeric",length=0) TMP2<-vector("numeric",length=0) TMP3<-vector("numeric",length=0) TMP4<-vector("numeric",length=0) TMP5<-vector("numeric",length=0) TMP6<-vector("numeric",length=0) TMP7<-vector("numeric",length=0) TMP8<-vector("numeric",length=0) for (lambda in lambdas){ IlWC<-INT - lambda*IWc IlWA<-INT - lambda*IWa IlWCinv <- solve(IlWC) IlWAinv <- solve(IlWA) tmp1<-vector("numeric",length=samples) tmp2<-vector("numeric",length=samples) tmp3<-vector("numeric",length=samples) tmp4<-vector("numeric",length=samples) tmp5<-vector("numeric",length=samples) tmp6<-vector("numeric",length=samples) tmp7<-vector("numeric",length=samples) tmp8<-vector("numeric",length=samples) i=1 while(i<=samples){ vit0 <- rnorm(NT,mean = 0, sd = 1) vit1 <- rnorm(NT,mean = 0, sd = 1) vit0nn <- (exp(rnorm(NT,0,1)) - exp(0.5)) /sqrt(exp(2)-exp(1)) vit1nn <- (exp(rnorm(NT,0,1)) - exp(0.5)) /sqrt(exp(2)-exp(1)) uit0 <- IrWCinv %*% vit0 uit1 <- vit1 uit0nn <- IrWCinv %*% vit0nn uit1nn <- vit1nn y0<- as.matrix(IlWCinv %*% (xtotb + uit0)) y1<- as.matrix(IlWAinv %*% (xtotb + uit1)) y0nn <- as.matrix(IlWCinv %*% (xtotb + uit0nn)) y1nn <- as.matrix(IlWAinv %*% (xtotb + uit1nn)) colnames(y0) <- "y0" colnames(y1) <- "y1" colnames(y0nn) <- "y0nn" colnames(y1nn) <- "y1nn" data1 <- data.frame(cbind(index, y0, y1, y0nn, y1nn, xtot)) H0_size <- spgm(y0 ~ x11 + x12, data = data1, index = NULL, listw = Wc, model = "within", lag = T) H1_size <- spgm(y0 ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H0_power <- spgm(y1 ~ x11 + x12, data = data1, index = NULL, listw = Wc, model = "within", lag = T) H1_power <- spgm(y1 ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error =F, method = "w2sls") H0_sizenn <- spgm(y0nn ~ x11 + x12, data = data1, index = NULL, listw = Wc, model = "within", lag = T) H1_sizenn <- spgm(y0nn ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H0_powernn <- spgm(y1nn ~ x11 + x12, data = data1, index = NULL, listw = Wc, model = "within", lag = T) H1_powernn <- spgm(y1nn ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error =F, method = "w2sls") pred_size<-predictors(H0_size, H1_size, y = y0, xreg = xtot, w0 = IWc, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_power<-predictors(H0_power, H1_power, y = y1, xreg = xtot, w0 = IWc, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_sizenn<-predictors(H0_sizenn, H1_sizenn, y = y0nn, xreg = xtot, w0 = IWc, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_powernn<-predictors(H0_powernn, H1_powernn, y = y1nn, xreg = xtot, w0 = IWc, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) test_size_p1<-try(kptest(y0,H0_size, pred = pred_size[,1], x0 = xtot, w0 = IWc, instruments = Hexp2, nobs = nobs, time = time)) test_size_p2<-try(kptest(y0, H0_size, pred = pred_size[,2], x0 = xtot, w0 = IWc, instruments = Hexp2, nobs = nobs, time = time)) test_power_p1<-try(kptest(y1,H0_power, pred=pred_power[,1], x0 = xtot, w0 = IWc, instruments = Hexp2, nobs = nobs, time = time)) test_power_p2<-try(kptest(y1,H0_power, pred=pred_power[,2], x0 = xtot, w0 = IWc, instruments = Hexp2, nobs = nobs, time = time)) test_size_p1nn<-try(kptest(y0nn,H0_sizenn, pred = pred_sizenn[,1], x0 = xtot, w0 = IWc, instruments = Hexp2, nobs = nobs, time = time)) test_size_p2nn<-try(kptest(y0nn, H0_sizenn, pred = pred_sizenn[,2], x0 = xtot, w0 = IWc, instruments = Hexp2, nobs = nobs, time = time)) test_power_p1nn<-try(kptest(y1nn,H0_powernn, pred=pred_powernn[,1], x0 = xtot, w0 = IWc, instruments = Hexp2, nobs = nobs, time = time)) test_power_p2nn<-try(kptest(y1nn,H0_powernn, pred=pred_powernn[,2], x0 = xtot, w0 = IWc, instruments = Hexp2, nobs = nobs, time = time)) tmp1[i]<-test_size_p1 tmp2[i]<-test_size_p2 tmp3[i]<-test_power_p1 tmp4[i]<-test_power_p2 tmp5[i]<-test_size_p1nn tmp6[i]<-test_size_p2nn tmp7[i]<-test_power_p1nn tmp8[i]<-test_power_p2nn i=i+1 if (i==samples) print(i) } TMP1<-c(TMP1, tmp1) TMP2<-c(TMP2, tmp2) TMP3<-c(TMP3, tmp3) TMP4<-c(TMP4, tmp4) TMP5<-c(TMP5, tmp5) TMP6<-c(TMP6, tmp6) TMP7<-c(TMP7, tmp7) TMP8<-c(TMP8, tmp8) } FIN1<-c(FIN1, TMP1) FIN2<-c(FIN2, TMP2) FIN3<-c(FIN3, TMP3) FIN4<-c(FIN4, TMP4) FIN5<-c(FIN5, TMP5) FIN6<-c(FIN6, TMP6) FIN7<-c(FIN7, TMP7) FIN8<-c(FIN8, TMP8) } STAT<- qchisq(0.05,1,lower.tail=FALSE) for (i in 1:length(FIN1)){ if(FIN1[i]>STAT) FIN1[i]=1 else FIN1[i]=0 if(FIN2[i]>STAT) FIN2[i]=1 else FIN2[i]=0 if(FIN3[i]>STAT) FIN3[i]=1 else FIN3[i]=0 if(FIN4[i]>STAT) FIN4[i]=1 else FIN4[i]=0 if(FIN5[i]>STAT) FIN5[i]=1 else FIN5[i]=0 if(FIN6[i]>STAT) FIN6[i]=1 else FIN6[i]=0 if(FIN7[i]>STAT) FIN7[i]=1 else FIN7[i]=0 if(FIN8[i]>STAT) FIN8[i]=1 else FIN8[i]=0 } org1<-matrix(FIN1,samples,length(lambdas)*length(rhos)) org2<-matrix(FIN2,samples,length(lambdas)*length(rhos)) org3<-matrix(FIN3,samples,length(lambdas)*length(rhos)) org4<-matrix(FIN4,samples,length(lambdas)*length(rhos)) org5<-matrix(FIN5,samples,length(lambdas)*length(rhos)) org6<-matrix(FIN6,samples,length(lambdas)*length(rhos)) org7<-matrix(FIN7,samples,length(lambdas)*length(rhos)) org8<-matrix(FIN8,samples,length(lambdas)*length(rhos)) FR1<-colSums(org1)/samples FR2<-colSums(org2)/samples FR3<-colSums(org3)/samples FR4<-colSums(org4)/samples FR5<-colSums(org5)/samples FR6<-colSums(org6)/samples FR7<-colSums(org7)/samples FR8<-colSums(org8)/samples table<-cbind(FR1,FR2, FR3, FR4, FR5,FR6, FR7, FR8) print(table) save(FIN1,FIN2,FIN3, FIN4, FIN5, FIN6, FIN7, FIN8, file = paste("exp2_",nobs,"_R2=",format(Rtilde, digits =2),sep="")) write.table(format(table, digits = 3, trim =TRUE),paste("tab2_",nobs,"_R2=",format(Rtilde, digits =2),".csv",sep=""), quote = FALSE, sep="&") ############################################### # Third Set of Experiments ############################################### FIN1<-vector("numeric",length=0) FIN2<-vector("numeric",length=0) FIN3<-vector("numeric",length=0) FIN4<-vector("numeric",length=0) FIN5<-vector("numeric",length=0) FIN6<-vector("numeric",length=0) FIN7<-vector("numeric",length=0) FIN8<-vector("numeric",length=0) for(rho in rhos){ IrWA<-INT - rho*IWa IrWAinv<-solve(IrWA) TMP1<-vector("numeric",length=0) TMP2<-vector("numeric",length=0) TMP3<-vector("numeric",length=0) TMP4<-vector("numeric",length=0) TMP5<-vector("numeric",length=0) TMP6<-vector("numeric",length=0) TMP7<-vector("numeric",length=0) TMP8<-vector("numeric",length=0) for (lambda in lambdas){ IlWA<-INT - lambda*IWa IlWAinv<-solve(IlWA) tmp1<-vector("numeric",length=samples) tmp2<-vector("numeric",length=samples) tmp3<-vector("numeric",length=samples) tmp4<-vector("numeric",length=samples) tmp5<-vector("numeric",length=samples) tmp6<-vector("numeric",length=samples) tmp7<-vector("numeric",length=samples) tmp8<-vector("numeric",length=samples) i=1 while(i<=samples){ vit0 <- rnorm(NT,mean = 0, sd = 1) vit1 <- rnorm(NT,mean = 0, sd = 1) vit0nn <- (exp(rnorm(NT,0,1)) - exp(0.5)) /sqrt(exp(2)-exp(1)) vit1nn <- (exp(rnorm(NT,0,1)) - exp(0.5)) /sqrt(exp(2)-exp(1)) uit0 <- IrWAinv %*% vit0 uit1 <- vit1 uit0nn <- IrWAinv %*% vit0nn uit1nn <- vit1nn y0<- as.matrix(IlWAinv %*% (xtotb + uit0)) y1<- as.matrix(IlWAinv %*% (xtot1b + uit1)) y0nn <- as.matrix(IlWAinv %*% (xtotb + uit0nn)) y1nn <- as.matrix(IlWAinv %*% (xtot1b + uit1nn)) colnames(y0) <- "y0" colnames(y1) <- "y1" colnames(y0nn) <- "y0nn" colnames(y1nn) <- "y1nn" data1 <- data.frame(cbind(index, y0, y1, y0nn, y1nn, xtot, xtot1)) H0_size <- spgm(y0 ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T) H1_size <- spgm(y0 ~ x21 + x22, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H0_power <- spgm(y1 ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T) H1_power <- spgm(y1 ~ x21 + x22, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error =F, method = "w2sls") H0_sizenn <- spgm(y0nn ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T) H1_sizenn <- spgm(y0nn ~ x21 + x22, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H0_powernn <- spgm(y1nn ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T) H1_powernn <- spgm(y1nn ~ x21 + x22, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error =F, method = "w2sls") pred_size<-predictors(H0_size, H1_size, y = y0, xreg = xtot1, w0 = IWa, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_power<-predictors(H0_power, H1_power, y = y1, xreg = xtot1, w0 = IWa, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_sizenn<-predictors(H0_sizenn, H1_sizenn, y = y0nn, xreg = xtot1, w0 = IWa, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_powernn<-predictors(H0_powernn, H1_powernn, y = y1nn, xreg = xtot1, w0 = IWa, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) test_size_p1<-try(kptest(y0,H0_size, pred = pred_size[,1], x0 = xtot, w0 = IWa, instruments = Hexp3, nobs = nobs, time = time)) test_size_p2<-try(kptest(y0, H0_size, pred = pred_size[,2], x0 = xtot, w0 = IWa, instruments = Hexp3, nobs = nobs, time = time)) test_power_p1<-try(kptest(y1,H0_power, pred=pred_power[,1], x0 = xtot, w0 = IWa, instruments = Hexp3, nobs = nobs, time = time)) test_power_p2<-try(kptest(y1,H0_power, pred=pred_power[,2], x0 = xtot, w0 = IWa, instruments = Hexp3, nobs = nobs, time = time)) test_size_p1nn<-try(kptest(y0nn,H0_sizenn, pred = pred_sizenn[,1], x0 = xtot, w0 = IWa, instruments = Hexp3, nobs = nobs, time = time)) test_size_p2nn<-try(kptest(y0nn, H0_sizenn, pred = pred_sizenn[,2], x0 = xtot, w0 = IWa, instruments = Hexp3, nobs = nobs, time = time)) test_power_p1nn<-try(kptest(y1nn,H0_powernn, pred=pred_powernn[,1], x0 = xtot, w0 = IWa, instruments = Hexp3, nobs = nobs, time = time)) test_power_p2nn<-try(kptest(y1nn,H0_powernn, pred=pred_powernn[,2], x0 = xtot, w0 = IWa, instruments = Hexp3, nobs = nobs, time = time)) tmp1[i]<-test_size_p1 tmp2[i]<-test_size_p2 tmp3[i]<-test_power_p1 tmp4[i]<-test_power_p2 tmp5[i]<-test_size_p1nn tmp6[i]<-test_size_p2nn tmp7[i]<-test_power_p1nn tmp8[i]<-test_power_p2nn i=i+1 if (i==samples) print(i) } TMP1<-c(TMP1, tmp1) TMP2<-c(TMP2, tmp2) TMP3<-c(TMP3, tmp3) TMP4<-c(TMP4, tmp4) TMP5<-c(TMP5, tmp5) TMP6<-c(TMP6, tmp6) TMP7<-c(TMP7, tmp7) TMP8<-c(TMP8, tmp8) } FIN1<-c(FIN1, TMP1) FIN2<-c(FIN2, TMP2) FIN3<-c(FIN3, TMP3) FIN4<-c(FIN4, TMP4) FIN5<-c(FIN5, TMP5) FIN6<-c(FIN6, TMP6) FIN7<-c(FIN7, TMP7) FIN8<-c(FIN8, TMP8) } STAT<- qchisq(0.05,1,lower.tail=FALSE) for (i in 1:length(FIN1)){ if(FIN1[i]>STAT) FIN1[i]=1 else FIN1[i]=0 if(FIN2[i]>STAT) FIN2[i]=1 else FIN2[i]=0 if(FIN3[i]>STAT) FIN3[i]=1 else FIN3[i]=0 if(FIN4[i]>STAT) FIN4[i]=1 else FIN4[i]=0 if(FIN5[i]>STAT) FIN5[i]=1 else FIN5[i]=0 if(FIN6[i]>STAT) FIN6[i]=1 else FIN6[i]=0 if(FIN7[i]>STAT) FIN7[i]=1 else FIN7[i]=0 if(FIN8[i]>STAT) FIN8[i]=1 else FIN8[i]=0 } org1<-matrix(FIN1,samples,length(lambdas)*length(rhos)) org2<-matrix(FIN2,samples,length(lambdas)*length(rhos)) org3<-matrix(FIN3,samples,length(lambdas)*length(rhos)) org4<-matrix(FIN4,samples,length(lambdas)*length(rhos)) org5<-matrix(FIN5,samples,length(lambdas)*length(rhos)) org6<-matrix(FIN6,samples,length(lambdas)*length(rhos)) org7<-matrix(FIN7,samples,length(lambdas)*length(rhos)) org8<-matrix(FIN8,samples,length(lambdas)*length(rhos)) FR1<-colSums(org1)/samples FR2<-colSums(org2)/samples FR3<-colSums(org3)/samples FR4<-colSums(org4)/samples FR5<-colSums(org5)/samples FR6<-colSums(org6)/samples FR7<-colSums(org7)/samples FR8<-colSums(org8)/samples table<-cbind(FR1,FR2, FR3, FR4, FR5,FR6, FR7, FR8) print(table) save(FIN1,FIN2,FIN3, FIN4, FIN5, FIN6, FIN7, FIN8, file = paste("exp3_",nobs,"_R2=",format(Rtilde, digits =2),sep="")) write.table(format(table, digits = 3, trim =TRUE),paste("tab3_",nobs,"_R2=",format(Rtilde, digits =2),".csv",sep=""), quote = FALSE, sep="&") ############################################### # Fourth Set of Experiments ############################################### FIN1<-vector("numeric",length=0) FIN2<-vector("numeric",length=0) FIN3<-vector("numeric",length=0) FIN4<-vector("numeric",length=0) FIN1bis<-vector("numeric",length=0) FIN2bis<-vector("numeric",length=0) FIN3bis<-vector("numeric",length=0) FIN4bis<-vector("numeric",length=0) FIN5<-vector("numeric",length=0) FIN6<-vector("numeric",length=0) FIN7<-vector("numeric",length=0) FIN8<-vector("numeric",length=0) FIN5bis<-vector("numeric",length=0) FIN6bis<-vector("numeric",length=0) FIN7bis<-vector("numeric",length=0) FIN8bis<-vector("numeric",length=0) for(rho in rhos){ IrWA<-INT - rho*IWa IrWC<-INT - rho*IWc IrWAinv<-solve(IrWA) IrWCinv<-solve(IrWC) TMP1<-vector("numeric",length=0) TMP2<-vector("numeric",length=0) TMP3<-vector("numeric",length=0) TMP4<-vector("numeric",length=0) TMP1bis<-vector("numeric",length=0) TMP2bis<-vector("numeric",length=0) TMP3bis<-vector("numeric",length=0) TMP4bis<-vector("numeric",length=0) TMP5<-vector("numeric",length=0) TMP6<-vector("numeric",length=0) TMP7<-vector("numeric",length=0) TMP8<-vector("numeric",length=0) TMP5bis<-vector("numeric",length=0) TMP6bis<-vector("numeric",length=0) TMP7bis<-vector("numeric",length=0) TMP8bis<-vector("numeric",length=0) for (lambda in lambdas){ IlWA<-INT - lambda*IWa IlWC<-INT - lambda*IWc IlWAinv<-solve(IlWA) IlWCinv<-solve(IlWC) tmp1<-vector("numeric",length=samples) tmp2<-vector("numeric",length=samples) tmp3<-vector("numeric",length=samples) tmp4<-vector("numeric",length=samples) tmp1bis<-vector("numeric",length=samples) tmp2bis<-vector("numeric",length=samples) tmp3bis<-vector("numeric",length=samples) tmp4bis<-vector("numeric",length=samples) tmp5<-vector("numeric",length=samples) tmp6<-vector("numeric",length=samples) tmp7<-vector("numeric",length=samples) tmp8<-vector("numeric",length=samples) tmp5bis<-vector("numeric",length=samples) tmp6bis<-vector("numeric",length=samples) tmp7bis<-vector("numeric",length=samples) tmp8bis<-vector("numeric",length=samples) i=1 while(i<=samples){ vit0 <- rnorm(NT,mean = 0, sd = 1) vit1 <- rnorm(NT,mean = 0, sd = 1) vit0nn <- (exp(rnorm(NT,0,1)) - exp(0.5)) /sqrt(exp(2)-exp(1)) vit1nn <- (exp(rnorm(NT,0,1)) - exp(0.5)) /sqrt(exp(2)-exp(1)) uit0<-IrWAinv %*% vit0 uit1<-IrWCinv %*% vit1 uit0nn<-IrWAinv %*% vit0nn uit1nn<-IrWCinv %*% vit1nn y0<- as.matrix(IlWAinv %*% (xtotb + uit0)) y1<- as.matrix(IlWCinv %*% (xtot1b + uit1)) y0nn <- as.matrix(IlWAinv %*% (xtotb + uit0nn)) y1nn <- as.matrix(IlWCinv %*% (xtot1b + uit1nn)) colnames(y0) <- "y0" colnames(y1) <- "y1" colnames(y0nn) <- "y0nn" colnames(y1nn) <- "y1nn" data1 <- data.frame(cbind(index, y0, y1, y0nn, y1nn, xtot)) H0_size <- spgm(y0 ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T) H1_size <- spgm(y0 ~ x21 + x22, data = data1, index = NULL, listw = Wc, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H1_sizebis <- spgm(y0 ~ x21 + x22, data = data1, index = NULL, listw = Wc, model = "within", lag = T, spatial.error = T) # # simply 2sls H0_power <- spgm(y1 ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T) H1_power <- spgm(y1 ~ x21 + x22, data = data1, index = NULL, listw = Wc, model = "within", lag = T, spatial.error =F, method = "w2sls") H1_powerbis <- spgm(y1 ~ x21 + x22, data = data1, index = NULL, listw = Wc, model = "within", lag = T, spatial.error =T) H0_sizenn <- spgm(y0nn ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T) H1_sizenn <- spgm(y0nn ~ x21 + x22, data = data1, index = NULL, listw = Wc, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H1_sizennbis <- spgm(y0nn ~ x21 + x22, data = data1, index = NULL, listw = Wc, model = "within", lag = T, spatial.error = T) H0_powernn <- spgm(y1nn ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T) H1_powernn <- spgm(y1nn ~ x21 + x22, data = data1, index = NULL, listw = Wc, model = "within", lag = T, spatial.error =F, method = "w2sls") H1_powernnbis <- spgm(y1nn ~ x21 + x22, data = data1, index = NULL, listw = Wc, model = "within", lag = T, spatial.error =T) pred_size<-predictors(H0_size, H1_size, y = y0, xreg = xtot1, w0 = IWa, w1 = IWc, nobs = nobs, time = time, IN = IN, IT = Itime) pred_power<-predictors(H0_power, H1_power, y = y1, xreg = xtot1, w0 = IWa, w1 = IWc, nobs = nobs, time = time, IN = IN, IT = Itime) pred_sizebis <-predictors(H0_size, H1_sizebis, y = y0, xreg = xtot1, w0 = IWa, w1 = IWc, nobs = nobs, time = time, IN = IN, IT = Itime) pred_powerbis <-predictors(H0_power, H1_powerbis, y = y1, xreg = xtot1, w0 = IWa, w1 = IWc, nobs = nobs, time = time, IN = IN, IT = Itime) pred_sizenn<-predictors(H0_sizenn, H1_sizenn, y = y0nn, xreg = xtot1, w0 = IWa, w1 = IWc, nobs = nobs, time = time, IN = IN, IT = Itime) pred_powernn<-predictors(H0_powernn, H1_powernn, y = y1nn, xreg = xtot1, w0 = IWa, w1 = IWc, nobs = nobs, time = time, IN = IN, IT = Itime) pred_sizennbis<-predictors(H0_sizenn, H1_sizennbis, y = y0nn, xreg = xtot1, w0 = IWa, w1 = IWc, nobs = nobs, time = time, IN = IN, IT = Itime) pred_powernnbis<-predictors(H0_powernn, H1_powernnbis, y = y1nn, xreg = xtot1, w0 = IWa, w1 = IWc, nobs = nobs, time = time, IN = IN, IT = Itime) test_size_p1<-try(kptest(y0,H0_size, pred = pred_size[,1], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_size_p2<-try(kptest(y0, H0_size, pred = pred_size[,2], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_size_p1bis<-try(kptest(y0,H0_size, pred = pred_sizebis[,1], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_size_p2bis<-try(kptest(y0, H0_size, pred = pred_sizebis[,2], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_power_p1<-try(kptest(y1,H0_power, pred=pred_power[,1], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_power_p2<-try(kptest(y1,H0_power, pred=pred_power[,2], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_power_p1bis<-try(kptest(y1,H0_power, pred=pred_powerbis[,1], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_power_p2bis<-try(kptest(y1,H0_power, pred=pred_powerbis[,2], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_size_p1nn<-try(kptest(y0nn,H0_sizenn, pred = pred_sizenn[,1], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_size_p2nn<-try(kptest(y0nn, H0_sizenn, pred = pred_sizenn[,2], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_size_p1nnbis<-try(kptest(y0nn,H0_sizenn, pred = pred_sizennbis[,1], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_size_p2nnbis<-try(kptest(y0nn, H0_sizenn, pred = pred_sizennbis[,2], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_power_p1nn<-try(kptest(y1nn,H0_powernn, pred=pred_powernn[,1], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_power_p2nn<-try(kptest(y1nn,H0_powernn, pred=pred_powernn[,2], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_power_p1nnbis<-try(kptest(y1nn,H0_powernn, pred=pred_powernnbis[,1], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) test_power_p2nnbis<-try(kptest(y1nn,H0_powernn, pred=pred_powernnbis[,2], x0 = xtot, w0 = IWa, instruments = Hexp4, nobs = nobs, time = time)) tmp1[i]<-test_size_p1 tmp2[i]<-test_size_p2 tmp3[i]<-test_power_p1 tmp4[i]<-test_power_p2 tmp1bis[i]<-test_size_p1bis tmp2bis[i]<-test_size_p2bis tmp3bis[i]<-test_power_p1bis tmp4bis[i]<-test_power_p2bis tmp5[i]<-test_size_p1nn tmp6[i]<-test_size_p2nn tmp7[i]<-test_power_p1nn tmp8[i]<-test_power_p2nn tmp5bis[i]<-test_size_p1nnbis tmp6bis[i]<-test_size_p2nnbis tmp7bis[i]<-test_power_p1nnbis tmp8bis[i]<-test_power_p2nnbis i=i+1 if (i==samples) print(i) } TMP1<-c(TMP1, tmp1) TMP2<-c(TMP2, tmp2) TMP3<-c(TMP3, tmp3) TMP4<-c(TMP4, tmp4) TMP1bis<-c(TMP1bis, tmp1bis) TMP2bis<-c(TMP2bis, tmp2bis) TMP3bis<-c(TMP3bis, tmp3bis) TMP4bis<-c(TMP4bis, tmp4bis) TMP5<-c(TMP5, tmp5) TMP6<-c(TMP6, tmp6) TMP7<-c(TMP7, tmp7) TMP8<-c(TMP8, tmp8) TMP5bis<-c(TMP5bis, tmp5bis) TMP6bis<-c(TMP6bis, tmp6bis) TMP7bis<-c(TMP7bis, tmp7bis) TMP8bis<-c(TMP8bis, tmp8bis) } FIN1<-c(FIN1, TMP1) FIN2<-c(FIN2, TMP2) FIN3<-c(FIN3, TMP3) FIN4<-c(FIN4, TMP4) FIN1bis<-c(FIN1bis, TMP1bis) FIN2bis<-c(FIN2bis, TMP2bis) FIN3bis<-c(FIN3bis, TMP3bis) FIN4bis<-c(FIN4bis, TMP4bis) FIN5<-c(FIN5, TMP5) FIN6<-c(FIN6, TMP6) FIN7<-c(FIN7, TMP7) FIN8<-c(FIN8, TMP8) FIN5bis<-c(FIN5bis, TMP5bis) FIN6bis<-c(FIN6bis, TMP6bis) FIN7bis<-c(FIN7bis, TMP7bis) FIN8bis<-c(FIN8bis, TMP8bis) } STAT<- qchisq(0.05,1,lower.tail=FALSE) for (i in 1:length(FIN1)){ if(FIN1[i]>STAT) FIN1[i]=1 else FIN1[i]=0 if(FIN1bis[i]>STAT) FIN1bis[i]=1 else FIN1bis[i]=0 if(FIN2[i]>STAT) FIN2[i]=1 else FIN2[i]=0 if(FIN2bis[i]>STAT) FIN2bis[i]=1 else FIN2bis[i]=0 if(FIN3[i]>STAT) FIN3[i]=1 else FIN3[i]=0 if(FIN3bis[i]>STAT) FIN3bis[i]=1 else FIN3bis[i]=0 if(FIN4[i]>STAT) FIN4[i]=1 else FIN4[i]=0 if(FIN4bis[i]>STAT) FIN4bis[i]=1 else FIN4bis[i]=0 if(FIN5[i]>STAT) FIN5[i]=1 else FIN5[i]=0 if(FIN5bis[i]>STAT) FIN5bis[i]=1 else FIN5bis[i]=0 if(FIN6[i]>STAT) FIN6[i]=1 else FIN6[i]=0 if(FIN6bis[i]>STAT) FIN6bis[i]=1 else FIN6bis[i]=0 if(FIN7[i]>STAT) FIN7[i]=1 else FIN7[i]=0 if(FIN7bis[i]>STAT) FIN7bis[i]=1 else FIN7bis[i]=0 if(FIN8[i]>STAT) FIN8[i]=1 else FIN8[i]=0 if(FIN8bis[i]>STAT) FIN8bis[i]=1 else FIN8bis[i]=0 } org1<-matrix(FIN1,samples,length(lambdas)*length(rhos)) org2<-matrix(FIN2,samples,length(lambdas)*length(rhos)) org3<-matrix(FIN3,samples,length(lambdas)*length(rhos)) org4<-matrix(FIN4,samples,length(lambdas)*length(rhos)) org1bis<-matrix(FIN1bis,samples,length(lambdas)*length(rhos)) org2bis<-matrix(FIN2bis,samples,length(lambdas)*length(rhos)) org3bis<-matrix(FIN3bis,samples,length(lambdas)*length(rhos)) org4bis<-matrix(FIN4bis,samples,length(lambdas)*length(rhos)) org5<-matrix(FIN5,samples,length(lambdas)*length(rhos)) org6<-matrix(FIN6,samples,length(lambdas)*length(rhos)) org7<-matrix(FIN7,samples,length(lambdas)*length(rhos)) org8<-matrix(FIN8,samples,length(lambdas)*length(rhos)) org5bis<-matrix(FIN5bis,samples,length(lambdas)*length(rhos)) org6bis<-matrix(FIN6bis,samples,length(lambdas)*length(rhos)) org7bis<-matrix(FIN7bis,samples,length(lambdas)*length(rhos)) org8bis<-matrix(FIN8bis,samples,length(lambdas)*length(rhos)) FR1<-colSums(org1)/samples FR2<-colSums(org2)/samples FR3<-colSums(org3)/samples FR4<-colSums(org4)/samples FR1bis<-colSums(org1bis)/samples FR2bis<-colSums(org2bis)/samples FR3bis<-colSums(org3bis)/samples FR4bis<-colSums(org4bis)/samples FR5<-colSums(org5)/samples FR6<-colSums(org6)/samples FR7<-colSums(org7)/samples FR8<-colSums(org8)/samples FR5bis<-colSums(org5bis)/samples FR6bis<-colSums(org6bis)/samples FR7bis<-colSums(org7bis)/samples FR8bis<-colSums(org8bis)/samples table<-cbind(FR1,FR1bis, FR2, FR2bis, FR3, FR3bis, FR4,FR4bis, FR5, FR5bis, FR6,FR6bis, FR7,FR7bis, FR8, FR8bis) print(table) save(FIN1,FIN1bis,FIN2,FIN2bis,FIN3,FIN3bis, FIN4,FIN4bis, FIN5,FIN5bis, FIN6,FIN6bis, FIN7,FIN7bis, FIN8,FIN8bis, file=paste("exp4_",nobs,"_R2=",format(Rtilde, digits =2),sep="") ) write.table(format(table, digits = 3, trim =TRUE),paste("tab4_",nobs,"_R2=",format(Rtilde, digits =2),".csv",sep=""), quote = FALSE, sep="&") ############################################## # Fifth Set of Experiments ############################################## FIN1<-vector("numeric",length=0) FIN2<-vector("numeric",length=0) FIN3<-vector("numeric",length=0) FIN4<-vector("numeric",length=0) FIN5<-vector("numeric",length=0) FIN6<-vector("numeric",length=0) FIN7<-vector("numeric",length=0) FIN8<-vector("numeric",length=0) for(rho in rhos){ IrWA<-INT - rho*IWa IrWC<-INT - rho*IWc IrWAinv<-solve(IrWA) IrWCinv<-solve(IrWC) TMP1<-vector("numeric",length=0) TMP2<-vector("numeric",length=0) TMP3<-vector("numeric",length=0) TMP4<-vector("numeric",length=0) TMP5<-vector("numeric",length=0) TMP6<-vector("numeric",length=0) TMP7<-vector("numeric",length=0) TMP8<-vector("numeric",length=0) for (lambda in lambdas){ IlWA<-INT - lambda*IWa IlWC<-INT - lambda*IWc IlWAinv<-solve(IlWA) IlWCinv<-solve(IlWC) tmp1<-vector("numeric",length=samples) tmp2<-vector("numeric",length=samples) tmp3<-vector("numeric",length=samples) tmp4<-vector("numeric",length=samples) tmp5<-vector("numeric",length=samples) tmp6<-vector("numeric",length=samples) tmp7<-vector("numeric",length=samples) tmp8<-vector("numeric",length=samples) i=1 while(i<=samples){ vit0 <- rnorm(NT,mean = 0, sd = 1) vit1 <- rnorm(NT,mean = 0, sd = 1) vit0nn <- (exp(rnorm(NT,0,1)) - exp(0.5)) /sqrt(exp(2)-exp(1)) vit1nn <- (exp(rnorm(NT,0,1)) - exp(0.5)) /sqrt(exp(2)-exp(1)) uit0 <- IrWCinv %*% vit0 uit1 <- vit1 uit0nn <- IrWCinv %*% vit0nn uit1nn <- vit1nn y0<- as.matrix(IlWCinv %*% (xtotb + uit0)) y1<- as.matrix(IlWAinv %*% (xtotb + uit1)) y0nn <- as.matrix(IlWCinv %*% (xtotb + uit0nn)) y1nn <- as.matrix(IlWAinv %*% (xtotb + uit1nn)) colnames(y0) <- "y0" colnames(y1) <- "y1" colnames(y0nn) <- "y0nn" colnames(y1nn) <- "y1nn" data1 <- data.frame(cbind(index, y0, y1, y0nn, y1nn, xtot)) H0_size <- spgm(y0 ~ x11 + x12, data = data1, index = NULL, listw = Wc, model = "within", lag = T) H1_size <- spgm(y0 ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H1_size2 <- spgm(y0 ~ x11 + x12, data = data1, index = NULL, listw = Wb, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H0_power <- spgm(y1 ~ x11 + x12, data = data1, index = NULL, listw = Wc, model = "within", lag = T) H1_power <- spgm(y1 ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error =F, method = "w2sls") H1_power2 <- spgm(y1 ~ x11 + x12, data = data1, index = NULL, listw = Wb, model = "within", lag = T, spatial.error =F, method = "w2sls") H0_sizenn <- spgm(y0nn ~ x11 + x12, data = data1, index = NULL, listw = Wc, model = "within", lag = T) H1_sizenn <- spgm(y0nn ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H1_sizenn2 <- spgm(y0nn ~ x11 + x12, data = data1, index = NULL, listw = Wb, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H0_powernn <- spgm(y1nn ~ x11 + x12, data = data1, index = NULL, listw = Wc, model = "within", lag = T) H1_powernn <- spgm(y1nn ~ x11 + x12, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error =F, method = "w2sls") H1_powernn2 <- spgm(y1nn ~ x11 + x12, data = data1, index = NULL, listw = Wb, model = "within", lag = T, spatial.error =F, method = "w2sls") pred_size<-predictors(H0_size, H1_size, y = y0, xreg = xtot, w0 = IWc, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_size2<-predictors(H0_size, H1_size, y = y0, xreg = xtot, w0 = IWc, w1 = IWb, nobs = nobs, time = time, IN = IN, IT = Itime) pred_power<-predictors(H0_power, H1_power, y = y1, xreg = xtot, w0 = IWc, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_power2<-predictors(H0_power, H1_power, y = y1, xreg = xtot, w0 = IWc, w1 = IWb, nobs = nobs, time = time, IN = IN, IT = Itime) pred_sizenn<-predictors(H0_sizenn, H1_sizenn, y = y0nn, xreg = xtot, w0 = IWc, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_sizenn2<-predictors(H0_sizenn, H1_sizenn, y = y0nn, xreg = xtot, w0 = IWc, w1 = IWb, nobs = nobs, time = time, IN = IN, IT = Itime) pred_powernn<-predictors(H0_powernn, H1_powernn, y = y1nn, xreg = xtot, w0 = IWc, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_powernn2<-predictors(H0_powernn, H1_powernn, y = y1nn, xreg = xtot, w0 = IWc, w1 = IWb, nobs = nobs, time = time, IN = IN, IT = Itime) test_size_p1<-try(kptest(y0,H0_size, pred = cbind(pred_size[,1], pred_size2[,1]), x0 = xtot, w0 = IWc, instruments = Hexp5, nobs = nobs, time = time)) test_size_p2<-try(kptest(y0, H0_size, pred = cbind(pred_size[,2], pred_size2[,2]), x0 = xtot, w0 = IWc, instruments = Hexp5, nobs = nobs, time = time)) test_power_p1<-try(kptest(y1,H0_power, pred=cbind(pred_power[,1],pred_power2[,1]), x0 = xtot, w0 = IWc, instruments = Hexp5, nobs = nobs, time = time)) test_power_p2<-try(kptest(y1,H0_power, pred=cbind(pred_power[,2],pred_power2[,2]), x0 = xtot, w0 = IWc, instruments = Hexp5, nobs = nobs, time = time)) test_size_p1nn<-try(kptest(y0nn,H0_sizenn, pred = cbind(pred_sizenn[,1],pred_sizenn2[,1]), x0 = xtot, w0 = IWc, instruments = Hexp5, nobs = nobs, time = time)) test_size_p2nn<-try(kptest(y0nn, H0_sizenn, pred = cbind(pred_sizenn[,2],pred_sizenn2[,2]), x0 = xtot, w0 = IWc, instruments = Hexp5, nobs = nobs, time = time)) test_power_p1nn<-try(kptest(y1nn,H0_powernn, pred = cbind(pred_powernn[,1],pred_powernn2[,1]), x0 = xtot, w0 = IWc, instruments = Hexp5, nobs = nobs, time = time)) test_power_p2nn<-try(kptest(y1nn,H0_powernn, pred = cbind(pred_powernn[,2],pred_powernn2[,2]), x0 = xtot, w0 = IWc, instruments = Hexp5, nobs = nobs, time = time)) tmp1[i]<-test_size_p1 tmp2[i]<-test_size_p2 tmp3[i]<-test_power_p1 tmp4[i]<-test_power_p2 tmp5[i]<-test_size_p1nn tmp6[i]<-test_size_p2nn tmp7[i]<-test_power_p1nn tmp8[i]<-test_power_p2nn i=i+1 if (i==samples) print(i) } TMP1<-c(TMP1, tmp1) TMP2<-c(TMP2, tmp2) TMP3<-c(TMP3, tmp3) TMP4<-c(TMP4, tmp4) TMP5<-c(TMP5, tmp5) TMP6<-c(TMP6, tmp6) TMP7<-c(TMP7, tmp7) TMP8<-c(TMP8, tmp8) } FIN1<-c(FIN1, TMP1) FIN2<-c(FIN2, TMP2) FIN3<-c(FIN3, TMP3) FIN4<-c(FIN4, TMP4) FIN5<-c(FIN5, TMP5) FIN6<-c(FIN6, TMP6) FIN7<-c(FIN7, TMP7) FIN8<-c(FIN8, TMP8) } STAT<- qchisq(0.05,2,lower.tail=FALSE) for (i in 1:length(FIN1)){ if(FIN1[i]>STAT) FIN1[i]=1 else FIN1[i]=0 if(FIN2[i]>STAT) FIN2[i]=1 else FIN2[i]=0 if(FIN3[i]>STAT) FIN3[i]=1 else FIN3[i]=0 if(FIN4[i]>STAT) FIN4[i]=1 else FIN4[i]=0 if(FIN5[i]>STAT) FIN5[i]=1 else FIN5[i]=0 if(FIN6[i]>STAT) FIN6[i]=1 else FIN6[i]=0 if(FIN7[i]>STAT) FIN7[i]=1 else FIN7[i]=0 if(FIN8[i]>STAT) FIN8[i]=1 else FIN8[i]=0 } org1<-matrix(FIN1,samples,length(lambdas)*length(rhos)) org2<-matrix(FIN2,samples,length(lambdas)*length(rhos)) org3<-matrix(FIN3,samples,length(lambdas)*length(rhos)) org4<-matrix(FIN4,samples,length(lambdas)*length(rhos)) org5<-matrix(FIN5,samples,length(lambdas)*length(rhos)) org6<-matrix(FIN6,samples,length(lambdas)*length(rhos)) org7<-matrix(FIN7,samples,length(lambdas)*length(rhos)) org8<-matrix(FIN8,samples,length(lambdas)*length(rhos)) FR1<-colSums(org1)/samples FR2<-colSums(org2)/samples FR3<-colSums(org3)/samples FR4<-colSums(org4)/samples FR5<-colSums(org5)/samples FR6<-colSums(org6)/samples FR7<-colSums(org7)/samples FR8<-colSums(org8)/samples table<-cbind(FR1,FR2, FR3, FR4, FR5,FR6, FR7, FR8) print(table) save(FIN1,FIN2,FIN3, FIN4, FIN5, FIN6, FIN7, FIN8, file = paste("exp5_",nobs,"_R2=",format(Rtilde, digits =2),sep="")) write.table(format(table, digits = 3, trim =TRUE),paste("tab5_",nobs,"_R2=",format(Rtilde, digits =2),".csv",sep=""), quote = FALSE, sep="&") ################################################ # Sixth Set of Experiments ############################################### FIN1<-vector("numeric",length=0) FIN2<-vector("numeric",length=0) FIN3<-vector("numeric",length=0) FIN4<-vector("numeric",length=0) FIN5<-vector("numeric",length=0) FIN6<-vector("numeric",length=0) FIN7<-vector("numeric",length=0) FIN8<-vector("numeric",length=0) for(rho in rhos){ IrWA<-INT - rho*IWa IrWC<-INT - rho*IWc IrWAinv<-solve(IrWA) IrWCinv<-solve(IrWC) TMP1<-vector("numeric",length=0) TMP2<-vector("numeric",length=0) TMP3<-vector("numeric",length=0) TMP4<-vector("numeric",length=0) TMP5<-vector("numeric",length=0) TMP6<-vector("numeric",length=0) TMP7<-vector("numeric",length=0) TMP8<-vector("numeric",length=0) for (lambda in lambdas){ IlWA<-INT - lambda*IWa IlWC<-INT - lambda*IWc IlWAinv<-solve(IlWA) IlWCinv<-solve(IlWC) tmp1<-vector("numeric",length=samples) tmp2<-vector("numeric",length=samples) tmp3<-vector("numeric",length=samples) tmp4<-vector("numeric",length=samples) tmp5<-vector("numeric",length=samples) tmp6<-vector("numeric",length=samples) tmp7<-vector("numeric",length=samples) tmp8<-vector("numeric",length=samples) i=1 while(i<=samples){ vit0 <- rnorm(NT,mean = 0, sd = 1) vit1 <- rnorm(NT,mean = 0, sd = 1) vit0nn <- (exp(rnorm(NT,0,1)) - exp(0.5)) /sqrt(exp(2)-exp(1)) vit1nn <- (exp(rnorm(NT,0,1)) - exp(0.5)) /sqrt(exp(2)-exp(1)) uit0 <- IrWCinv %*% vit0 uit1 <- vit1 uit0nn <- IrWCinv %*% vit0nn uit1nn <- vit1nn y0<- as.matrix(IlWCinv %*% (xtotb + uit0)) y1<- as.matrix(IlWAinv %*% (xtot1b + uit1)) y0nn <- as.matrix(IlWCinv %*% (xtotb + uit0nn)) y1nn <- as.matrix(IlWAinv %*% (xtot1b + uit1nn)) colnames(y0) <- "y0" colnames(y1) <- "y1" colnames(y0nn) <- "y0nn" colnames(y1nn) <- "y1nn" data1 <- data.frame(cbind(index, y0, y1, y0nn, y1nn, xtot, xtot1)) H0_size <- spgm(y0 ~ x11 + x12, data = data1, index = NULL, listw = Wc, model = "within", lag = T) H1_size <- spgm(y0 ~ x21 + x22, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H1_size2 <- spgm(y0 ~ x21 + x22, data = data1, index = NULL, listw = Wb, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H0_power <- spgm(y1 ~ x11 + x12, data = data1, index = NULL, listw = Wc, model = "within", lag = T) H1_power <- spgm(y1 ~ x21 + x22, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error =F, method = "w2sls") H1_power2 <- spgm(y1 ~ x21 + x22, data = data1, index = NULL, listw = Wb, model = "within", lag = T, spatial.error =F, method = "w2sls") H0_sizenn <- spgm(y0nn ~ x11 + x12, data = data1, index = NULL, listw = Wc, model = "within", lag = T) H1_sizenn <- spgm(y0nn ~ x21 + x22, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H1_sizenn2 <- spgm(y0nn ~ x21 + x22, data = data1, index = NULL, listw = Wb, model = "within", lag = T, spatial.error = F, method = "w2sls") # # simply 2sls H0_powernn <- spgm(y1nn ~ x11 + x12, data = data1, index = NULL, listw = Wc, model = "within", lag = T) H1_powernn <- spgm(y1nn ~ x21 + x22, data = data1, index = NULL, listw = Wa, model = "within", lag = T, spatial.error =F, method = "w2sls") H1_powernn2 <- spgm(y1nn ~ x21 + x22, data = data1, index = NULL, listw = Wb, model = "within", lag = T, spatial.error =F, method = "w2sls") pred_size<-predictors(H0_size, H1_size, y = y0, xreg = xtot1, w0 = IWc, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_size2<-predictors(H0_size, H1_size, y = y0, xreg = xtot1, w0 = IWc, w1 = IWb, nobs = nobs, time = time, IN = IN, IT = Itime) pred_power<-predictors(H0_power, H1_power, y = y1, xreg = xtot1, w0 = IWc, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_power2<-predictors(H0_power, H1_power, y = y1, xreg = xtot1, w0 = IWc, w1 = IWb, nobs = nobs, time = time, IN = IN, IT = Itime) pred_sizenn<-predictors(H0_sizenn, H1_sizenn, y = y0nn, xreg = xtot1, w0 = IWc, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_sizenn2<-predictors(H0_sizenn, H1_sizenn, y = y0nn, xreg = xtot1, w0 = IWc, w1 = IWb, nobs = nobs, time = time, IN = IN, IT = Itime) pred_powernn<-predictors(H0_powernn, H1_powernn, y = y1nn, xreg = xtot1, w0 = IWc, w1 = IWa, nobs = nobs, time = time, IN = IN, IT = Itime) pred_powernn2<-predictors(H0_powernn, H1_powernn, y = y1nn, xreg = xtot1, w0 = IWc, w1 = IWb, nobs = nobs, time = time, IN = IN, IT = Itime) test_size_p1<-try(kptest(y0,H0_size, pred = cbind(pred_size[,1], pred_size2[,1]), x0 = xtot, w0 = IWc, instruments = Hexp6, nobs = nobs, time = time)) test_size_p2<-try(kptest(y0, H0_size, pred = cbind(pred_size[,2], pred_size2[,2]), x0 = xtot, w0 = IWc, instruments = Hexp6, nobs = nobs, time = time)) test_power_p1<-try(kptest(y1,H0_power, pred=cbind(pred_power[,1],pred_power2[,1]), x0 = xtot, w0 = IWc, instruments = Hexp6, nobs = nobs, time = time)) test_power_p2<-try(kptest(y1,H0_power, pred=cbind(pred_power[,2],pred_power2[,2]), x0 = xtot, w0 = IWc, instruments = Hexp6, nobs = nobs, time = time)) test_size_p1nn<-try(kptest(y0nn,H0_sizenn, pred = cbind(pred_sizenn[,1],pred_sizenn2[,1]), x0 = xtot, w0 = IWc, instruments = Hexp6, nobs = nobs, time = time)) test_size_p2nn<-try(kptest(y0nn, H0_sizenn, pred = cbind(pred_sizenn[,2],pred_sizenn2[,2]), x0 = xtot, w0 = IWc, instruments = Hexp6, nobs = nobs, time = time)) test_power_p1nn<-try(kptest(y1nn,H0_powernn, pred = cbind(pred_powernn[,1],pred_powernn2[,1]), x0 = xtot, w0 = IWc, instruments = Hexp6, nobs = nobs, time = time)) test_power_p2nn<-try(kptest(y1nn,H0_powernn, pred = cbind(pred_powernn[,2],pred_powernn2[,2]), x0 = xtot, w0 = IWc, instruments = Hexp6, nobs = nobs, time = time)) tmp1[i]<-test_size_p1 tmp2[i]<-test_size_p2 tmp3[i]<-test_power_p1 tmp4[i]<-test_power_p2 tmp5[i]<-test_size_p1nn tmp6[i]<-test_size_p2nn tmp7[i]<-test_power_p1nn tmp8[i]<-test_power_p2nn i=i+1 if (i==samples) print(i) } TMP1<-c(TMP1, tmp1) TMP2<-c(TMP2, tmp2) TMP3<-c(TMP3, tmp3) TMP4<-c(TMP4, tmp4) TMP5<-c(TMP5, tmp5) TMP6<-c(TMP6, tmp6) TMP7<-c(TMP7, tmp7) TMP8<-c(TMP8, tmp8) } FIN1<-c(FIN1, TMP1) FIN2<-c(FIN2, TMP2) FIN3<-c(FIN3, TMP3) FIN4<-c(FIN4, TMP4) FIN5<-c(FIN5, TMP5) FIN6<-c(FIN6, TMP6) FIN7<-c(FIN7, TMP7) FIN8<-c(FIN8, TMP8) } STAT<- qchisq(0.05,2,lower.tail=FALSE) for (i in 1:length(FIN1)){ if(FIN1[i]>STAT) FIN1[i]=1 else FIN1[i]=0 if(FIN2[i]>STAT) FIN2[i]=1 else FIN2[i]=0 if(FIN3[i]>STAT) FIN3[i]=1 else FIN3[i]=0 if(FIN4[i]>STAT) FIN4[i]=1 else FIN4[i]=0 if(FIN5[i]>STAT) FIN5[i]=1 else FIN5[i]=0 if(FIN6[i]>STAT) FIN6[i]=1 else FIN6[i]=0 if(FIN7[i]>STAT) FIN7[i]=1 else FIN7[i]=0 if(FIN8[i]>STAT) FIN8[i]=1 else FIN8[i]=0 } org1<-matrix(FIN1,samples,length(lambdas)*length(rhos)) org2<-matrix(FIN2,samples,length(lambdas)*length(rhos)) org3<-matrix(FIN3,samples,length(lambdas)*length(rhos)) org4<-matrix(FIN4,samples,length(lambdas)*length(rhos)) org5<-matrix(FIN5,samples,length(lambdas)*length(rhos)) org6<-matrix(FIN6,samples,length(lambdas)*length(rhos)) org7<-matrix(FIN7,samples,length(lambdas)*length(rhos)) org8<-matrix(FIN8,samples,length(lambdas)*length(rhos)) FR1<-colSums(org1)/samples FR2<-colSums(org2)/samples FR3<-colSums(org3)/samples FR4<-colSums(org4)/samples FR5<-colSums(org5)/samples FR6<-colSums(org6)/samples FR7<-colSums(org7)/samples FR8<-colSums(org8)/samples table<-cbind(FR1,FR2, FR3, FR4, FR5,FR6, FR7, FR8) print(table) save(FIN1,FIN2,FIN3, FIN4, FIN5, FIN6, FIN7, FIN8, file = paste("exp6_",nobs,"_R2=",format(Rtilde, digits =2),sep="")) write.table(format(table, digits = 3, trim =TRUE),paste("tab6_",nobs,"_R2=",format(Rtilde, digits =2),".csv",sep=""), quote = FALSE, sep="&") }