-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathipr_generated_from_fake_series.R
76 lines (57 loc) · 2.83 KB
/
ipr_generated_from_fake_series.R
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
rm(list=ls())
require(ggplot2)
require(ggthemes)
require(statnet)
require(coda)
require(igraph)
require(foreach)
require(pracma)
require(rgexf)
require(fields)
require(data.table)
require(igraph)
require(doParallel)
numTimeSeries_Space = round(linspace(2,50,n=10))
lengthOfTimeSeries_Space = round(linspace(20,1000,n=10))
test_frame = data.table(expand.grid(numTimeSeries_Space, lengthOfTimeSeries_Space))
colnames(test_frame) = c("NUM_SERIES", "LENGTH_SERIES")
test_frame = test_frame[LENGTH_SERIES > NUM_SERIES]
numSamples = 1E4
iprSimple = function(inputMatrix) { return (colSums(eigen(inputMatrix)$vectors^4)) }
num_cores = detectCores() - 1
cl <- makeCluster(num_cores)
registerDoParallel(cl)
clusterExport(cl, c("numSamples","test_frame", "iprSimple"))
allIPRks <- foreach(rc = 1:nrow(test_frame), .packages=c("data.table")) %dopar%
{
p = test_frame[rc, NUM_SERIES]
n = test_frame[rc, LENGTH_SERIES]
S = diag(p)
draws = rWishart(numSamples, n, S)/n
iprks = t(apply(draws, 3, iprSimple))
return(list("NUM_SERIES"=p,
"LENGTH_SERIES"=n,
"IPR_k"=iprks))
}
IPRKSims <- lapply( lapply (allIPRks, "[[", "IPR_k"),
function(x) { structure(x, class = c("IPRKSim","matrix"))} )
iprkFrame = data.table("NUM_SERIES" = unlist(lapply(allIPRks, "[[", "NUM_SERIES")),
"LENGTH_SERIES" = unlist(lapply(allIPRks, "[[", "LENGTH_SERIES")),
"IPRKs" = IPRKSims
)
plot(density(unlist(iprkFrame[NUM_SERIES == 2 & LENGTH_SERIES == 129, IPRKs])))
plot(density(unlist(iprkFrame[NUM_SERIES == 7 & LENGTH_SERIES == 129, IPRKs])))
plot(density(unlist(iprkFrame[NUM_SERIES == 23 & LENGTH_SERIES == 129, IPRKs])))
plot(density(unlist(iprkFrame[NUM_SERIES == 50 & LENGTH_SERIES == 129, IPRKs])))
plot(density(unlist(iprkFrame[NUM_SERIES == 50 & LENGTH_SERIES == 1000, IPRKs])))
plot(density(unlist(iprkFrame[NUM_SERIES == 7 & LENGTH_SERIES == 1000, IPRKs])))
plot(density(unlist(iprkFrame[NUM_SERIES == 2 & LENGTH_SERIES == 1000, IPRKs])))
plot(density(unlist(iprkFrame[NUM_SERIES == 50 & LENGTH_SERIES == 456, IPRKs])))
plot(density(unlist(iprkFrame[NUM_SERIES == 7 & LENGTH_SERIES == 456, IPRKs])))
plot(density(unlist(iprkFrame[NUM_SERIES == 2 & LENGTH_SERIES == 456, IPRKs])))
# just individual IPR_k where k=1,10
plot(density(unlist(iprkFrame[NUM_SERIES == 18 & LENGTH_SERIES == 1000, IPRKs][[1]][,1] )))
plot(density(unlist(iprkFrame[NUM_SERIES == 18 & LENGTH_SERIES == 1000, IPRKs][[1]][,10] )))
plot(density(unlist(iprkFrame[NUM_SERIES == 18 & LENGTH_SERIES == 1000, IPRKs][[1]][,18] )))
plot(density(unlist(iprkFrame[NUM_SERIES == 2 & LENGTH_SERIES == 1000, IPRKs][[1]][,1] )))
plot(density(unlist(iprkFrame[NUM_SERIES == 2 & LENGTH_SERIES == 1000, IPRKs][[1]][,2] )))