In [1]:

    

include("ilmtools.jl")


    

    












    












    
Out[1]:





create_loglikelihood (generic function with 1 method)

In [2]:

    

pop_db1 = create_pop_db_univariates(100, Uniform(0,10), Uniform(0,10))


    

    
Out[2]:




ind_id
x
y
1
1
7.321472426772846
0.706273939122668
2
2
9.231480924022065
9.846036948667118
3
3
1.3830527445491336
3.8892082096289737
4
4
5.942387136217954
7.254722392731283
5
5
4.680869648483172
8.33989052973654
6
6
3.69939239500102
4.598482547418767
7
7
0.832671120854751
9.729997231932558
8
8
3.132337833945338
9.136637241904534
9
9
4.137903525533771
1.631541477712981
10
10
6.711171512847769
8.365907203913466
11
11
1.3423059065884613
6.484382473805206
12
12
1.8674494910056016
6.430064041751066
13
13
4.059220614923607
8.587607500935835
14
14
9.72157663683841
7.493942118752061
15
15
2.074569844127463
4.7177304848878565
16
16
1.562840404146233
0.7574889241288107
17
17
3.8881081617965774
0.36735118330215144
18
18
2.8620543892371964
6.752863189195772
19
19
5.09660137661045
8.657344078629945
20
20
6.810707012857757
8.711124368957407
21
21
7.842049808670229
6.91163257557416
22
22
7.742471257051302
6.3302570984383095
23
23
1.3400618338672299
2.3208863054534645
24
24
8.436071015336815
5.691154670125229
25
25
6.98616233235259
3.0670702324665355
26
26
4.524275955641879
5.133520268081531
27
27
4.475728289153887
1.9133377748825242
28
28
3.4948089846645547
8.607843597666742
29
29
1.9544921138956162
0.5594929107777058
30
30
1.1673802095511077
7.442061912761442
&vellip
&vellip
&vellip
&vellip

In [3]:

    

evdb = infect_recover_loop(pop_db1, "discrete", "SIR", 2, 6)


    

    
Out[3]:





edb(100x3 DataFrame
| Row | ind_id | s | i    |
|-----|--------|---|------|
| 1   | 1      | 0 | 8.0  |
| 2   | 2      | 0 | 16.0 |
| 3   | 3      | 0 | 5.0  |
| 4   | 4      | 0 | 9.0  |
| 5   | 5      | 0 | 9.0  |
| 6   | 6      | 0 | 6.0  |
| 7   | 7      | 0 | 12.0 |
| 8   | 8      | 0 | 11.0 |
| 9   | 9      | 0 | 5.0  |
| 10  | 10     | 0 | 11.0 |
| 11  | 11     | 0 | 9.0  |
⋮
| 89  | 89     | 0 | 13.0 |
| 90  | 90     | 0 | 9.0  |
| 91  | 91     | 0 | 4.0  |
| 92  | 92     | 0 | 14.0 |
| 93  | 93     | 0 | 15.0 |
| 94  | 94     | 0 | 9.0  |
| 95  | 95     | 0 | 11.0 |
| 96  | 96     | 0 | 8.0  |
| 97  | 97     | 0 | 11.0 |
| 98  | 98     | 0 | 5.0  |
| 99  | 99     | 0 | 4.0  |
| 100 | 100    | 0 | 8.0  |,[1.0,2.0,2.0,3.0,3.0,3.0,4.0,4.0,4.0,4.0  …  13.0,13.0,13.0,14.0,14.0,14.0,14.0,15.0,15.0,16.0],"discrete")

In [4]:

    

ll = create_loglikelihood(pop_db1, evdb)


    

    
Out[4]:





dSI_ll (generic function with 1 method)

In [5]:

    

ll((2, 6))
@time ll((2, 6))


    

    




elapsed time: 0.000992869 seconds (376984 bytes allocated)






    
Out[5]:





-73.74265521426908

In [6]:

    

n = 20000
burnin = 10000
sim = Chains(n, 2, names = ["alpha", "beta"])
theta = AMMVariate([2.0, 6.0])
SigmaF = cholfact(eye(2))
@time for i in 1:n
    amm!(theta, SigmaF, ll, adapt = (i <= burnin))
    sim[i,:,1] = theta[1:2]
end
Mamba.describe(sim)


    

    




elapsed time: 54.885969224 seconds (14565503944 bytes allocated, 27.72% gc time)
Iterations = 1:20000
Thinning interval = 1
Chains = 1
Samples per chain = 20000

Empirical Posterior Estimates:
3x6 Array{Any,2}:
 ""        "Mean"   "SD"      "Naive SE"   "MCSE"         "ESS"
 "alpha"  1.81213  0.374386  0.00264731   0.00688161  7693.86  
 "beta"   6.43645  0.514608  0.00363883   0.0101574   7164.86  

Quantiles:
3x6 Array{Any,2}:
 ""        "2.5%"   "25.0%"   "50.0%"   "75.0%"   "97.5%"
 "alpha"  1.17396  1.54992   1.7803    2.03924   2.635   
 "beta"   5.49502  6.06808   6.4071    6.77212   7.52396 

In [7]:

    

plot(x=burnin:n, y=sim.value[burnin:n, 1, 1], Geom.line)


    

    
Out[7]:

In [8]:

    

plot(x=burnin:n, y=sim.value[burnin:n, 2, 1], Geom.line)


    

    
Out[8]: