Bricolage.Rmd 1.47 KB
 Eric Coissac committed Jan 24, 2019 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 --- title: "R Notebook" output: html_document: df_print: paged --- {r} Xs = test2 nperm=100  {r} ix = c(4,3,6,5,8,7) iy = 1 r = corls(Xs,permutations = how(nperm = nperm)) rx = r[ix,ix] ry = r[iy,ix] rxi = r[ix[(ix %% 2)>0],ix[(ix %% 2)>0]] ryi = r[iy,ix[(ix %% 2)>0]]  {r} r2 = (solve(rx) %*% ry) * ry r2  {r} sum(r2)  {r} r2 = (solve(rxi) %*% ryi) * ryi r2  {r} sum(r2)  {r} v = varls(Xs,permutations = how(nperm = nperm)) vt = v[c(iy,ix[(ix %% 2)>0]),c(iy,ix[(ix %% 2)>0])] * (nrow(Xs)-1) vr = v[c(iy+1,ix[(ix %% 2)==0]),c(iy,ix[(ix %% 2)>0])] * (nrow(Xs)-1) vt vr  {r} #diag(vr)=0 vp = vt - vr vp  {r} d = sqrt(diag(vp)) n = d %o% d r = vp / n r  {r} rx = r[ix[(ix %% 2)==0]/2,ix[(ix %% 2)==0]/2] ry = r[iy,ix[(ix %% 2)==0]/2] r2 = (solve(rx) %*% ry) * ry r2  {r} sum(r2)  {r} rp = solve(r) rp = rp / sqrt(diag(rp) %o% diag(rp)) rp rp^2  {r} library(ProcMod) r2s = seq(from=0.1, to = 0.9, by = 0.1) n.sim=20 r2.obs = matrix(0,nrow = n.sim,ncol = length(r2s)) d1 = simulate_matrix(10,2000) for (j in seq_along(r2s)) { print(c("r2",as.character(r2s[j]))) for (i in seq_len(n.sim)) { print(c("round",as.character(i))) d2 = simulate_correlation(10,3000,d1,r2s[j]) v = varls2(d1,d2,permutations = 10) s = sqrt(diag(v)) r2.obs[i,j] = (v / outer(s,s))[1,2]^2 } } boxplot(r2.obs) cbind(R2_theo=r2s, R2_obs=colMeans(r2.obs), R2_sd=apply(r2.obs,MARGIN = 2,FUN = sd))