Commit c2b5a844 authored by Jorge Visca's avatar Jorge Visca

R_e computation

parent 8babd937
......@@ -5,8 +5,6 @@ local lib
local mathutil = require('lib.mathutil')
local people
local G = require('lib.graph')
G.load('in/nodes.csv', 'in/edges.csv')
......@@ -18,18 +16,19 @@ local function w_log (t, ...)
--conf.f_log:flush()
end
M.people = nil
-- Creates a new person, assigning processess
M.create_exposed = function ( zone, exposed_from )
local t = lib.esim.now
local new_exposed = {
id=#people+1,
id=#M.people+1,
exposed_from = exposed_from,
exposed_to = {},
in_zone = zone,
disease = lib.disease.new_disease()
}
people[new_exposed.id] = new_exposed
M.people[new_exposed.id] = new_exposed
if type(exposed_from)=='table' then
exposed_from.exposed_to[#exposed_from.exposed_to+1] = new_exposed
end
......@@ -85,7 +84,7 @@ M.initialize = function ( c )
end
mobility.initialize( conf )
people = {}
M.people = {}
if not conf.f_log then w_log = function() end end
end
......
......@@ -23,10 +23,16 @@ stats filename nooutput
# "out/states1.log" using 1:3 w steps, \
# "out/states1.log" using 1:4 w steps, \
plot for [i=3:STATS_columns] filename using 1:i w l
plot for [i=5:STATS_columns] filename using 1:i w l
set out "contacttracingR0.eps"
set ylabel 'R_e'
unset yrange
set yrange [0:]
plot filename using 1:2:3 w e t '', \
filename using 1:2 w l t '', \
unset output
local M = {}
local mathutil = require('lib.mathutil')
M.run = function ( conf )
local lib = {}
lib.esim = require('lib.esim')
......@@ -16,6 +18,8 @@ M.run = function ( conf )
--conf.f_log:flush()
end
local R_e_A = {}
--os.execute('rm out/*.log out/*.txt')
for i = 1, conf.NRUNS do
......@@ -54,17 +58,38 @@ M.run = function ( conf )
until sim_status
w_log(t_now, ' finished ',sim_status,' R: ', lib.disease.state_count.R)
print (i, 'finished at t: '..t_now..' events: '
..event_count..' ('..sim_status ..')')
if conf.f_states then conf.f_states:close() end
if conf.f_log then conf.f_log:close() end
sim_status_count[sim_status] = (sim_status_count[sim_status] or 0) + 1
print (i, 'finished at t: '..t_now..' events: '
..event_count..' ('..sim_status ..')')
do
local a = {}
for id, agent in pairs(lib.agents.people) do
if agent.disease.state == 'R' then
a[#a+1] = #agent.exposed_to
end
end
if #a>0 then
local R_e = mathutil.avgstd(a)
R_e_A[#R_e_A+1] = R_e
end
end
end
return sim_status_count
local R_e_avg, R_e_std = mathutil.avgstd(R_e_A)
local sim_result = {
status_count = sim_status_count,
R_e_avg = R_e_avg,
R_e_std = R_e_std,
}
return sim_result
end
return M
\ No newline at end of file
......@@ -23,6 +23,7 @@ local Efirst = {}
local Zhit = {}
local Ecount = {}
local ExtintionR = {}
local ReffA = {}
local FinishedT = {}
......@@ -31,13 +32,20 @@ local function parse (f)
local Es = {}
local Ef = {}
local status, R
local Infcount = {}
for l in f:lines() do
local t, event, s = l:match('^(%S+)%s(%S+)%s?(.*)$')
t = assert(tonumber(t))
if event=='new_exposed' then
--230.6 new_exposed id: 75 mobility: static from: 9 zone: 7 CCZ07
local zid, zone = s:match('zone: (%d+) (.*)$')
local id, iid, zid, zone = s:match('id: (%d+) mobility: %S+ from: (%S+) zone: (%d+) (.*)$')
id = tonumber(id)
Infcount[id] = 0
iid = tonumber(iid)
if iid then -- skipt infector 'none'
Infcount[iid] = Infcount[iid] + 1
end
local zonelabel = zone --zid
Es[zonelabel] = (Es[zonelabel] or 0) + 1
elseif event=='firstzone_exposed' then
......@@ -80,6 +88,15 @@ local function parse (f)
ZhitS[zone] = (ZhitS[zone] or 0) + 1
end
end
do
local IcA = {}
for id, count in pairs(Infcount) do
IcA[#IcA+1] = count
end
local Reff = mathutil.avgstd(IcA)
ReffA[#ReffA+1] = Reff
end
end
......@@ -271,6 +288,8 @@ end
parse_files()
print('R_e avg std:', mathutil.avgstd(ReffA))
process_finished()
process_Ecount()
process_Efirst()
......
......@@ -6,7 +6,7 @@ local epitrace = require('epitrace')
local conf = {
TMAX = math.huge, --12*30 -- Simulation time limit
NRUNS = 10000, -- Number of runs
NRUNS = 1000, -- Number of runs
out_path = 'out/', -- nil to disable logging
startExposed = {
......@@ -64,14 +64,17 @@ conf.scenario_process = function (conf)
end
end
local sim_status_count = epitrace.run(conf)
-- run simulation
local sim_result = epitrace.run(conf)
-- Print run totals
for status, count in pairs(sim_status_count) do
print('Re avg,std:', sim_result.R_e_avg, sim_result.R_e_std)
for status, count in pairs(sim_result.status_count) do
print('runs', status, count, count/conf.NRUNS)
end
-- Launch plotting (on background)
print('plotting launched in background')
os.execute('gnuplot states.plot &')
os.execute('luajit parse_logs.lua '..conf.NRUNS..' &')
......@@ -65,22 +65,25 @@ conf.scenario_process = function (conf)
end
local frunlog = io.open('out/contacttracing.log', 'w')
frunlog:write('tracers runout outbreak extinction hitstop\n')
frunlog:write('tracers Re ReStd runout outbreak extinction hitstop\n')
-- run simulations
for ntracers = 0, 20, 2 do
conf.contacttracers = ntracers
conf.contacttracers = ntracers -- set sim. scenario
local s = epitrace.run(conf)
local s = epitrace.run(conf) -- run
frunlog:write(
ntracers, ' ',
s.runout/conf.NRUNS, ' ',
s.outbreak/conf.NRUNS, ' ',
s.extinction/conf.NRUNS, ' ',
s.hit_stop/conf.NRUNS, '\n'
s.R_e_avg, ' ',
s.R_e_std, ' ',
(s.status_count.runout or 0)/conf.NRUNS, ' ',
(s.status_count.outbreak or 0)/conf.NRUNS, ' ',
(s.status_count.extinction or 0)/conf.NRUNS, ' ',
(s.status_count.hit_stop or 0)/conf.NRUNS, '\n'
)
end
frunlog:close()
os.execute('gnuplot contacttracing.plot &')
os.execute('gnuplot contacttracing.plot')
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment