Introduction to openMalariaUtilities
Source:vignettes/openmalariautilities.Rmd
openmalariautilities.Rmd
This document is still work in progress
OpenMalaria is a powerful tool for simulating and studying the epidemiology of malaria. This R package allows the use of OpenMalaria from R.
Prerequisites
Make sure that SQLite
and OpenMalaria
are
both installed and available on the system. Please note: Even though the
package should work on Windows, it is developed and tested under
Linux.
To install the package, run
devtools::install_github("SwissTPH/r-openMalariaUtilities")
Introduction
Within this tutorial, we will create and run example scenarios. Each scenario needs to be defined in a XML file which serves as an input for OpenMalaria (OM). OpenMalariaUtilities (OMU) contains many functions which help to generate these input files. In this tutorial we will re-create the example from OpenMalaria itself.
XML file generation
A large part of OMU is dedicated at the generation of valid OM XML files. Before we start with the generation of these files, we will explain the folder layout of your projects.
In the folder hierarchy, your working or project directory is the top level. Within this folder, you can have one or more folders which correspond to the names of experiments. The simulations of these experiments are controlled via one or multiple R scripts, placed in the project directory. Please note: Most of the directories can be modified via function arguments, but for the sake of simplicity we using the defaults here.
MyProject
├── experiment_1
├── experiment_2
├── ...
├── script_1.R
├── script_2.R
└── ...
Each of the experiment folders will have a similar layout and content.
MyProject/experiment_1
├── cache/
├── logs/
├── outputs/
├── scenarios/
├── autoRegressionParameters.csv
├── densities.csv
├── base.xml
└── scenario_44.xsd
You do not need to worry about the creation of the folders and files. OMU will create necessary folders and download required files for you.
The cache
folder contains R objects which are used by
OMU to store reusable information. This will speed up certain
operations. Log files are stored in the logs
folder (which
will contain sub-folders for each step in the process) and OpenMalaria’s
raw output will be outputs
. The XML files for each
simulation scenario are stored in the scenarios
folder.
Starting the project
In your R script, you will probably start with
library(openMalariaUtilities)
. Thus, the beginning could
look like
library(openMalariaUtilities)
## Basic skeleton
baseList <- list(
## Mandatory
expName = "example",
## Mandatory
OMVersion = 44L,
## Optional
## analysisNo = 1L,
## Mandatory
demography = list(),
monitoring = list(),
interventions = list(),
healthSystem = list(),
entomology = list(),
## These are optional for OM
## parasiteGenetics = list(),
## pharmacology = list(),
## diagnostics = list(),
model = list()
)
For XML file generation, OMU uses a list as input (here:
baseList
). This list has a very similar
structure as the XML file required by OM but is easier to edit from
R compared to XML. The above snippet corresponds to the minimal skeleton
of this input list. expName
is the name of your experiment
and OMVersion
the version of OpenMalaria. Currently, only
version 44 is supported. Also note, that the OM version needs to be an
integer.
Each of the six nested lists demography
,
monitoring
, interventions
,
healthSystem
, entomology
and
model
need to be filled for OM.
Demography
There are two ways to add the data to the demography
entry:
- Add it manually,
e.g.
baseList[["demography"]] <- list(name = "Ifakara", ...)
- Use the provided
defineDemography()
function
We will use the second option in this tutorial. You can find various
defineXY
functions which are each designed to add a
specific entry to the list. However, you always have the option to write
parts or change specific entries manually (See also
extractList()
).
demo_data <- data.frame(
poppercent = c(
3.474714994, 12.76004028, 14.52151394, 12.75565434, 10.83632374,
8.393312454, 7.001421452, 5.800587654, 5.102136612, 4.182561874,
3.339409351, 2.986112356, 2.555766582, 2.332763433, 1.77400255,
1.008525491, 0.74167341, 0.271863401, 0.161614642
),
upperbound = c(
1, 5, 10, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60, 65, 70, 75, 80, 85, 90
)
)
baseList <- defineDemography(
baseList,
name = "Ifakara",
popSize = 4000L,
maximumAgeYrs = 90,
lowerbound = 0,
poppercent = demo_data$poppercent,
upperbound = demo_data$upperbound
)
str(baseList, max.level = 2)
#> List of 8
#> $ expName : chr "example"
#> $ OMVersion : int 44
#> $ demography :List of 4
#> ..$ name : chr "Ifakara"
#> ..$ popSize : int 4000
#> ..$ maximumAgeYrs: num 90
#> ..$ ageGroup :List of 20
#> $ monitoring : list()
#> $ interventions: list()
#> $ healthSystem : list()
#> $ entomology : list()
#> $ model : list()
In the above snippet, we started by creating a data frame which
contains the demographic data. This in turn was used as input for
defineDemography()
among other parameters. The function
will take care that the information is added at the correct position and
also check, that the correct data types are used.
Monitoring
The monitoring
section controls the output of
OpenMalaria, which is generated during the conducted surveys.
baseList[["monitoring"]] <- list(
name = "Quarterly Surveys",
## Mandatory, different from OM schema
startDate = "2000-01-01",
continuous = monitoringContinuousGen(
period = 1,
options = list(
name = c("input EIR", "simulated EIR", "human infectiousness", "N_v0",
"immunity h", "immunity Y", "new infections",
"num transmitting humans", "ITN coverage", "GVI coverage", "alpha",
"P_B", "P_C*P_D"),
value = c("true", "true", "true", "true", "true", "true", "true", "true",
"true", "true", "true", "false", "false")
)
),
SurveyOptions = monitoringSurveyOptionsGen(
options = list(
name = c("nHost", "nPatent", "nUncomp", "nSevere", "nDirDeaths",
"inputEIR", "simulatedEIR"),
value = c("true", "true", "true", "true", "false", "true", "true")
)
)
)
We started the monitoring definition by giving it a name and defining a start date. Here, this is mandatory because it allows us to use dates as input for OpenMalaria. It defines the starting point of the survey period.
For the options, two generator functions are used which take the name of the option and its value as input.
baseList[["monitoring"]] <- c(
baseList[["monitoring"]],
list(
surveys = monitoringSurveyTimesGen(
detectionLimit = 40, startDate = "2000-01-01", endDate = "2020-01-01",
interval = "quarterly", simStart = "1980-01-01"
),
ageGroup = surveyAgeGroupsGen(
lowerbound = 0,
upperbounds = c(5, 90)
)
)
)
str(baseList[["monitoring"]], max.level = 1)
#> List of 6
#> $ name : chr "Quarterly Surveys"
#> $ startDate : chr "2000-01-01"
#> $ continuous :List of 14
#> $ SurveyOptions:List of 7
#> $ surveys :List of 86
#> $ ageGroup :List of 3
Similar to the age groups in the demography section,
surveyAgeGroupsGen()
adds the desired age groups for the
surveys.
The definition of the survey times itself is done above. Please note,
that you really should use monitoringSurveyTimesGen()
as
this function will make sure your desired survey points are valid (and
adjust them if necessary) and will store them for later use in the
cache.
We can inspect the cache to see if our survey times have been generated correctly.
listCache()
#>
#> Name: mon_ageGroups Timestamp: 2024-09-06 12:52:36
#> Class: list
#> Value:
#> $lowerbound
#> [1] 0
#>
#> $upperbounds
#> [1] 5 90
#>
#>
#> --------------------------------------------------------------------------------
#>
#> Name: simStart Timestamp: 2024-09-06 12:52:36
#> Class: character
#> Value:
#> [1] "1980-01-01"
#>
#> --------------------------------------------------------------------------------
#>
#> Name: surveyTimes Timestamp: 2024-09-06 12:52:36
#> Class: data.table data.frame
#> Value:
#> number date day timestep
#> <int> <Date> <num> <num>
#> 1: 1 1999-12-27 7295 1459
#> 2: 2 2000-03-27 7385 1477
#> 3: 3 2000-06-30 7480 1496
#> 4: 4 2000-09-28 7570 1514
#> 5: 5 2000-12-27 7660 1532
#> 6: 6 2001-03-27 7750 1550
#> 7: 7 2001-06-30 7845 1569
#> 8: 8 2001-09-28 7935 1587
#> 9: 9 2001-12-27 8025 1605
#> 10: 10 2002-03-27 8115 1623
#> 11: 11 2002-06-30 8210 1642
#> 12: 12 2002-09-28 8300 1660
#> 13: 13 2002-12-27 8390 1678
#> 14: 14 2003-03-27 8480 1696
#> 15: 15 2003-06-30 8575 1715
#> 16: 16 2003-09-28 8665 1733
#> 17: 17 2003-12-27 8755 1751
#> 18: 18 2004-03-27 8845 1769
#> 19: 19 2004-06-30 8940 1788
#> 20: 20 2004-09-28 9030 1806
#> 21: 21 2004-12-27 9120 1824
#> 22: 22 2005-03-27 9210 1842
#> 23: 23 2005-06-30 9305 1861
#> 24: 24 2005-09-28 9395 1879
#> 25: 25 2005-12-27 9485 1897
#> 26: 26 2006-03-27 9575 1915
#> 27: 27 2006-06-30 9670 1934
#> 28: 28 2006-09-28 9760 1952
#> 29: 29 2006-12-27 9850 1970
#> 30: 30 2007-03-27 9940 1988
#> 31: 31 2007-06-30 10035 2007
#> 32: 32 2007-09-28 10125 2025
#> 33: 33 2007-12-27 10215 2043
#> 34: 34 2008-03-27 10305 2061
#> 35: 35 2008-06-30 10400 2080
#> 36: 36 2008-09-28 10490 2098
#> 37: 37 2008-12-27 10580 2116
#> 38: 38 2009-03-27 10670 2134
#> 39: 39 2009-06-30 10765 2153
#> 40: 40 2009-09-28 10855 2171
#> 41: 41 2009-12-27 10945 2189
#> 42: 42 2010-03-27 11035 2207
#> 43: 43 2010-06-30 11130 2226
#> 44: 44 2010-09-28 11220 2244
#> 45: 45 2010-12-27 11310 2262
#> 46: 46 2011-03-27 11400 2280
#> 47: 47 2011-06-30 11495 2299
#> 48: 48 2011-09-28 11585 2317
#> 49: 49 2011-12-27 11675 2335
#> 50: 50 2012-03-27 11765 2353
#> 51: 51 2012-06-30 11860 2372
#> 52: 52 2012-09-28 11950 2390
#> 53: 53 2012-12-27 12040 2408
#> 54: 54 2013-03-27 12130 2426
#> 55: 55 2013-06-30 12225 2445
#> 56: 56 2013-09-28 12315 2463
#> 57: 57 2013-12-27 12405 2481
#> 58: 58 2014-03-27 12495 2499
#> 59: 59 2014-06-30 12590 2518
#> 60: 60 2014-09-28 12680 2536
#> 61: 61 2014-12-27 12770 2554
#> 62: 62 2015-03-27 12860 2572
#> 63: 63 2015-06-30 12955 2591
#> 64: 64 2015-09-28 13045 2609
#> 65: 65 2015-12-27 13135 2627
#> 66: 66 2016-03-27 13225 2645
#> 67: 67 2016-06-30 13320 2664
#> 68: 68 2016-09-28 13410 2682
#> 69: 69 2016-12-27 13500 2700
#> 70: 70 2017-03-27 13590 2718
#> 71: 71 2017-06-30 13685 2737
#> 72: 72 2017-09-28 13775 2755
#> 73: 73 2017-12-27 13865 2773
#> 74: 74 2018-03-27 13955 2791
#> 75: 75 2018-06-30 14050 2810
#> 76: 76 2018-09-28 14140 2828
#> 77: 77 2018-12-27 14230 2846
#> 78: 78 2019-03-27 14320 2864
#> 79: 79 2019-06-30 14415 2883
#> 80: 80 2019-09-28 14505 2901
#> 81: 81 2019-12-27 14595 2919
#> 82: 82 2020-03-27 14685 2937
#> 83: 83 2020-06-30 14780 2956
#> 84: 84 2020-09-28 14870 2974
#> 85: 85 2020-12-27 14960 2992
#> number date day timestep
#>
#> --------------------------------------------------------------------------------
We can notice two things:
- We are expecting 81 surveys, however we have 85
- The dates have been adjusted
OMU adds one more year (here: 4 surveys) to the survey period. This is done in order to ensure that all deployments show their effect. The dates are adjusted so they match OM’s internal time steps. For more details on this, see here.
Furthermore, the number of days or time steps is counted from the
simulation start simStart
, which is here the 1st January
1980. It is recommended to start the simulation 20 - 50 years before the
first survey.
Interventions
In this section we will define actual measures which will be used throughout the simulation.
baseList[["interventions"]] <- list(
name = "test",
human = list(
component = list(
id = "GVI",
name = "DDT test",
GVI = list(
decay = list(
L = "0.5",
"function" = "exponential"
),
anophelesParams = list(
mosquito = "gambiae_ss", propActive = "1",
deterrency = list(value = "0.56"),
preprandialKillingEffect = list(value = "0"),
postprandialKillingEffect = list(value = "0.24")
)
)
),
deployment = list(
name = "DDT test",
component = list(id = "GVI"),
timed = list(
deploy = list(coverage = "0.9", time = "10y")
)
)
)
)
Health system
baseList[["healthSystem"]] <- list(
ImmediateOutcomes = list(
name = "Tanzania ACT",
drugRegimen = list(
firstLine = "ACT",
inpatient = "QN",
secondLine = "ACT"
),
initialACR = list(
ACT = list(value = 0.85),
QN = list(value = 0.998),
selfTreatment = list(value = 0.63)
),
compliance = list(
ACT = list(value = 0.9),
selfTreatment = list(value = 0.85)
),
nonCompliersEffective = list(
ACT = list(value = 0),
selfTreatment = list(value = 0)
),
treatmentActions = list(
ACT = list(
name = "clear blood-stage infections",
clearInfections = list(
stage = "blood",
timesteps = "1"
)
),
QN = list(
name = "clear blood-stage infections",
clearInfections = list(
stage = "blood",
timesteps = "1"
)
)
),
pSeekOfficialCareUncomplicated1 = list(value = 0.04),
pSelfTreatUncomplicated = list(value = 0.01),
pSeekOfficialCareUncomplicated2 = list(value = 0.04),
pSeekOfficialCareSevere = list(value = 0.48)
),
CFR = list(
group = list(lowerbound = 0, value = 0.09189),
group = list(lowerbound = 0.25, value = 0.0810811),
group = list(lowerbound = 0.75, value = 0.0648649),
group = list(lowerbound = 1.5, value = 0.0689189),
group = list(lowerbound = 2.5, value = 0.0675676),
group = list(lowerbound = 3.5, value = 0.0297297),
group = list(lowerbound = 4.5, value = 0.0459459),
group = list(lowerbound = 7.5, value = 0.0945946),
group = list(lowerbound = 12.5, value = 0.1243243),
group = list(lowerbound = 15, value = 0.1378378)
),
## Mandatory
pSequelaeInpatient = list(
interpolation = "none",
group = list(lowerbound = "0.0", value = 0.0132),
group = list(lowerbound = "5.0", value = 0.005)
)
)
Entomology
baseList[["entomology"]] <- list(
mode = "dynamic", name = "Namawala",
vector = list(
anopheles = list(
mosquito = "gambiae_ss", propInfected = 0.078, propInfectious = "0.021",
seasonality = list(
annualEIR = "16", input = "EIR",
fourierSeries = list(
EIRRotateAngle = "0",
coeffic = list(
a = "0.8968", b = "2.678"
),
coeffic = list(a = "-0.4551", b = "2.599")
)
),
mosq = list(
minInfectedThreshold = "0.001",
mosqRestDuration = list(value = "3"),
extrinsicIncubationPeriod = list(value = "11"),
mosqLaidEggsSameDayProportion = list(value = "0.313"),
mosqSeekingDuration = list(value = "0.33"),
mosqSurvivalFeedingCycleProbability = list(value = "0.623"),
availability = list(distr = "const"),
mosqProbBiting = list(mean = "0.95", variance = "0"),
mosqProbFindRestSite = list(mean = "0.95", variance = "0"),
mosqProbResting = list(mean = "0.99", variance = "0"),
mosqProbOvipositing = list(value = "0.88"),
mosqHumanBloodIndex = list(value = "0.939")
),
nonHumanHosts = list(
name = "unprotectedAnimals",
mosqRelativeEntoAvailability = list(value = "1.0"),
mosqProbBiting = list(value = "0.95"),
mosqProbFindRestSite = list(value = "0.95"),
mosqProbResting = list(value = "0.99")
)
),
nonHumanHosts = list(name = "unprotectedAnimals", number = "1.0")
)
)
Model
baseList[["model"]] <- list(
ModelOptions = list(
option = list(name = "LOGNORMAL_MASS_ACTION", value = "true"),
option = list(name = "NO_PRE_ERYTHROCYTIC", value = "true"),
option = list(name = "INNATE_MAX_DENS", value = "false"),
option = list(name = "INDIRECT_MORTALITY_FIX", value = "false"),
option = list(name = "NON_MALARIA_FEVERS", value = "true")
),
clinical = list(
healthSystemMemory = "6",
NonMalariaFevers = list(
incidence = list(
group = list(lowerbound = "0", value = "0.322769924518357"),
group = list(lowerbound = "0", value = "0.308520194304172"),
group = list(lowerbound = "1", value = "0.279441774808493"),
group = list(lowerbound = "2", value = "0.250431781111273"),
group = list(lowerbound = "3", value = "0.223285859756841"),
group = list(lowerbound = "4", value = "0.199298352451799"),
group = list(lowerbound = "5", value = "0.179376872365614"),
group = list(lowerbound = "6", value = "0.163623659390782"),
group = list(lowerbound = "7", value = "0.152227726923469"),
group = list(lowerbound = "8", value = "0.145022785567758"),
group = list(lowerbound = "9", value = "0.141493087461765"),
group = list(lowerbound = "10", value = "0.140473293219353"),
group = list(lowerbound = "11", value = "0.141109775159515"),
group = list(lowerbound = "12", value = "0.142644475217328"),
group = list(lowerbound = "13", value = "0.144335079395766"),
group = list(lowerbound = "14", value = "0.145964032924869"),
group = list(lowerbound = "15", value = "0.147708915135714"),
group = list(lowerbound = "16", value = "0.149731543445568"),
group = list(lowerbound = "17", value = "0.151887428568276"),
group = list(lowerbound = "18", value = "0.154060663485195"),
group = list(lowerbound = "19", value = "0.156179169710494"),
group = list(lowerbound = "20", value = "0.158135015380583"),
group = list(lowerbound = "21", value = "0.159704766482219"),
group = list(lowerbound = "22", value = "0.160807788387655"),
group = list(lowerbound = "23", value = "0.161427976448279"),
group = list(lowerbound = "24", value = "0.161620429119137"),
group = list(lowerbound = "25", value = "0.16144021875986"),
group = list(lowerbound = "26", value = "0.160943264630612"),
group = list(lowerbound = "27", value = "0.160217573697398"),
group = list(lowerbound = "28", value = "0.159422614374451"),
group = list(lowerbound = "29", value = "0.158542519631641"),
group = list(lowerbound = "30", value = "0.157501217628248"),
group = list(lowerbound = "31", value = "0.156175160594841"),
group = list(lowerbound = "32", value = "0.154402302191411"),
group = list(lowerbound = "33", value = "0.152102040636481"),
group = list(lowerbound = "34", value = "0.14921450014676"),
group = list(lowerbound = "35", value = "0.145714433541659"),
group = list(lowerbound = "36", value = "0.141800502067518"),
group = list(lowerbound = "37", value = "0.137916853907569"),
group = list(lowerbound = "38", value = "0.134503529382102"),
group = list(lowerbound = "39", value = "0.131746276580642"),
group = list(lowerbound = "40", value = "0.12969902537497"),
group = list(lowerbound = "41", value = "0.128398077347679"),
group = list(lowerbound = "42", value = "0.127864136551891"),
group = list(lowerbound = "43", value = "0.12804497197004"),
group = list(lowerbound = "44", value = "0.128894055047661"),
group = list(lowerbound = "45", value = "0.130350838992718"),
group = list(lowerbound = "46", value = "0.132286605622701"),
group = list(lowerbound = "47", value = "0.134599921072495"),
group = list(lowerbound = "48", value = "0.137212726976988"),
group = list(lowerbound = "49", value = "0.140035253913284"),
group = list(lowerbound = "50", value = "0.142934573453621"),
group = list(lowerbound = "51", value = "0.145830221511879"),
group = list(lowerbound = "52", value = "0.148674810561069"),
group = list(lowerbound = "53", value = "0.151497963594518"),
group = list(lowerbound = "54", value = "0.15438856687865"),
group = list(lowerbound = "55", value = "0.157403790093505"),
group = list(lowerbound = "56", value = "0.16059513222516"),
group = list(lowerbound = "57", value = "0.16402433342886"),
group = list(lowerbound = "58", value = "0.16770481415944"),
group = list(lowerbound = "59", value = "0.171626873047865"),
group = list(lowerbound = "60", value = "0.175748327054247"),
group = list(lowerbound = "61", value = "0.180030857856799"),
group = list(lowerbound = "62", value = "0.184411365583771"),
group = list(lowerbound = "63", value = "0.188816421789366"),
group = list(lowerbound = "64", value = "0.19316997803338"),
group = list(lowerbound = "65", value = "0.197435603275487"),
group = list(lowerbound = "66", value = "0.201578808813379"),
group = list(lowerbound = "67", value = "0.205556806881398"),
group = list(lowerbound = "68", value = "0.209307183457343"),
group = list(lowerbound = "69", value = "0.212783260344084"),
group = list(lowerbound = "70", value = "0.215944154621391"),
group = list(lowerbound = "71", value = "0.218749275266548"),
group = list(lowerbound = "72", value = "0.221187990639016"),
group = list(lowerbound = "73", value = "0.223361260399378"),
group = list(lowerbound = "74", value = "0.225363436789592"),
group = list(lowerbound = "75", value = "0.227254280093211"),
group = list(lowerbound = "76", value = "0.229084576349576"),
group = list(lowerbound = "77", value = "0.230891971097789"),
group = list(lowerbound = "78", value = "0.232690225166173"),
group = list(lowerbound = "79", value = "0.234484973338876"),
group = list(lowerbound = "80", value = "0.236276361586796"),
group = list(lowerbound = "81", value = "0.238064394629696"),
group = list(lowerbound = "82", value = "0.239849077182917"),
group = list(lowerbound = "83", value = "0.241630413957381"),
group = list(lowerbound = "84", value = "0.243408409659591"),
group = list(lowerbound = "85", value = "0.245183068991633"),
group = list(lowerbound = "86", value = "0.246954396651183"),
group = list(lowerbound = "87", value = "0.248722397331501"),
group = list(lowerbound = "88", value = "0.250487075721441"),
group = list(lowerbound = "89", value = "0.252248436505447"),
group = list(lowerbound = "90", value = "0.253127874257909")
)
)
),
human = list(
availabilityToMosquitoes = list(
group = list(lowerbound = "0", value = "0.7076"),
group = list(lowerbound = "0", value = "0.8538"),
group = list(lowerbound = "5", value = "1.0"),
group = list(lowerbound = "5", value = "1.0")
)
),
parameters = list(
interval = "5", iseed = "1", latentp = "3",
parameter = list(include = "false", name = "'-ln(1-Sinf)'", number = "1", value = "0.050736"),
parameter = list(include = "false", name = "Estar", number = "2", value = "0.03247"),
parameter = list(include = "false", name = "Simm", number = "3", value = "0.138161050830301"),
parameter = list(include = "false", name = "Xstar_p", number = "4", value = "1514.385853233699891"),
parameter = list(include = "false", name = "gamma_p", number = "5", value = "2.03692533424484"),
parameter = list(include = "false", name = "sigma2i", number = "6", value = "10.173598698525799"),
parameter = list(include = "false", name = "CumulativeYstar", number = "7", value = "35158523.31132510304451"),
parameter = list(include = "false", name = "CumulativeHstar", number = "8", value = "97.334652723897705"),
parameter = list(include = "false", name = "'-ln(1-alpha_m)'", number = "9", value = "2.33031045876193"),
parameter = list(include = "false", name = "decay_m", number = "10", value = "2.53106547375805"),
parameter = list(include = "false", name = "sigma2_0", number = "11", value = "0.655747311168152"),
parameter = list(include = "false", name = "Xstar_v", number = "12", value = "0.916181104713054"),
parameter = list(include = "false", name = "Ystar2", number = "13", value = "6502.26335600001039"),
parameter = list(include = "false", name = "alpha", number = "14", value = "142601.912520000012591"),
parameter = list(include = "false", name = "Density bias (non Garki)", number = "15", value = "0.177378570987455"),
parameter = list(include = "false", name = " sigma2 ", number = "16", value = "0.05"),
parameter = list(include = "false", name = "log oddsr CF community", number = "17", value = "0.736202"),
parameter = list(include = "false", name = "Indirect risk cofactor", number = "18", value = "0.018777338"),
parameter = list(include = "false", name = "Non-malaria infant mortality", number = "19", value = "49.539046599999999"),
parameter = list(include = "false", name = "Density bias (Garki)", number = "20", value = "4.79610772546704"),
parameter = list(include = "false", name = "Severe Malaria Threshhold", number = "21", value = "784455.599999999976717"),
parameter = list(include = "false", name = "Immunity Penalty", number = "22", value = "1"),
parameter = list(include = "false", name = "Immune effector decay", number = "23", value = "0"),
parameter = list(include = "false", name = "comorbidity intercept", number = "24", value = "0.0968"),
parameter = list(include = "false", name = "Ystar half life", number = "25", value = "0.275437402"),
parameter = list(include = "false", name = "Ystar1", number = "26", value = "0.596539864"),
parameter = list(include = "false", name = "Asexual immunity decay", number = "27", value = "0"),
parameter = list(include = "false", name = "Ystar0", number = "28", value = "296.302437899999973"),
parameter = list(include = "false", name = "Idete multiplier", number = "29", value = "2.797523626"),
parameter = list(include = "false", name = "critical age for comorbidity", number = "30", value = "0.117383")
)
)
Prepare to run the simulations
## Setup dirs
setupDirs(experimentName = "exp_test", replace = TRUE)
## Copy necessary OpenMalaria files.
setupOM()
#> [1] TRUE
setupDirs()
creates the experiment folder given as an
argument in the project directory. This also includes all the
sub-folders which were mentioned in the introduction .
We can verify this by checking the layout of the experiment directory.
list.dirs("exp_test")
#> [1] "exp_test" "exp_test/cache"
#> [3] "exp_test/logs" "exp_test/logs/results"
#> [5] "exp_test/logs/scenarios" "exp_test/logs/simulation"
#> [7] "exp_test/outputs" "exp_test/scenarios"
setupOM()
will make sure that all required files for
OpenMalaria to run are in place. This includes the correct
.xsd
schema file (according to the specified OM version) as
well as further supplementary files.
list.files("exp_test")
#> [1] "autoRegressionParameters.csv" "cache"
#> [3] "densities.csv" "logs"
#> [5] "outputs" "scenario_44.xsd"
#> [7] "scenarios"
Create XML file
As the name implies, createBaseXml()
will translate our
input list to a XML document and write it to disk.
## Create base XML file
createBaseXml(baseList, replace = TRUE)
#> Warning in createBaseXml(baseList, replace = TRUE): startDate 2000-01-01 is not
#> the same as simStart 1980-01-01. This can cause survey times to not appear!
readLines(file("exp_test/exp_test_base.xml", "r"), n = 10)
#> [1] "<?xml version=\"1.0\" encoding=\"UTF-8\"?>"
#> [2] "<om:scenario xmlns:om=\"http://openmalaria.org/schema/scenario_44\" xmlns:xsi=\"http://www.w3.org/2001/XMLSchema-instance\" name=\"exp_test\" schemaVersion=\"44\" xsi:schemaLocation=\"http://openmalaria.org/schema/scenario_44 scenario_44.xsd\">"
#> [3] " <demography name=\"Ifakara\" popSize=\"4000\" maximumAgeYrs=\"90\">"
#> [4] " <ageGroup lowerbound=\"0\">"
#> [5] " <group poppercent=\"3.474714994\" upperbound=\"1\"/>"
#> [6] " <group poppercent=\"12.76004028\" upperbound=\"5\"/>"
#> [7] " <group poppercent=\"14.52151394\" upperbound=\"10\"/>"
#> [8] " <group poppercent=\"12.75565434\" upperbound=\"15\"/>"
#> [9] " <group poppercent=\"10.83632374\" upperbound=\"20\"/>"
#> [10] " <group poppercent=\"8.393312454\" upperbound=\"25\"/>"
Run simulations
In order to run the simulation, we can simply call
runSimulations()
.
## Run simulations
runSimulations()
#> [1] "Running scenario [1/1]"
By default, the output in the R console is fairly minimal. More
detailed logs can be found in the logs/simulation
folder.
list.files("exp_test/logs/simulation")
#> [1] "exp_test_base.log"
OpenMalaria generates two output
files which are both placed in the outputs
folder:
- The
*_cts.txt
file
read.table("exp_test/outputs/exp_test_base_cts.txt", header = TRUE, sep = "\t", skip = 1, nrows = 10)
#> timestep input.EIR simulated.EIR human.infectiousness N_v0.gambiae_ss.
#> 1 0 0.00995606 0.00995049 0.0179202 30147.4
#> 2 1 0.01958060 0.15390600 0.0175667 39181.0
#> 3 2 0.03850920 0.24658900 0.0171114 46368.5
#> 4 3 0.07460740 0.37191900 0.0167741 49838.5
#> 5 4 0.14030700 0.52390000 0.0164523 48641.0
#> 6 5 0.25256000 0.68747000 0.0165751 43197.7
#> 7 6 0.42953400 0.84148200 0.0188351 35064.0
#> 8 7 0.68226100 0.96738200 0.0231314 26185.7
#> 9 8 1.00214000 1.07333000 0.0296558 18146.4
#> 10 9 1.35047000 1.18202000 0.0377866 11790.1
#> immunity.h immunity.Y new.infections num.transmitting.humans ITN.coverage
#> 1 256.846 11899700 0 2613 0
#> 2 256.817 11900500 400 2590 0
#> 3 256.856 11895200 692 2572 0
#> 4 256.934 11896300 1054 2549 0
#> 5 257.109 11897600 1469 2526 0
#> 6 257.299 11899200 1795 2507 0
#> 7 257.791 11901800 2250 2636 0
#> 8 258.367 11903800 2643 2798 0
#> 9 259.013 11904600 2987 3019 0
#> 10 259.681 11909900 3256 3280 0
#> GVI.coverage alpha_i.gambiae_ss.
#> 1 0 0.000206347
#> 2 0 0.000206342
#> 3 0 0.000206349
#> 4 0 0.000206332
#> 5 0 0.000206344
#> 6 0 0.000206344
#> 7 0 0.000206347
#> 8 0 0.000206338
#> 9 0 0.000206330
#> 10 0 0.000206329
- The
*_out.txt
file
read.table("exp_test/outputs/exp_test_base_out.txt", header = FALSE, sep = "\t", nrows = 10)
#> V1 V2 V3 V4
#> 1 1 1 0 6.07000e+02
#> 2 1 2 0 3.39300e+03
#> 3 1 1 3 1.37000e+02
#> 4 1 2 3 1.15500e+03
#> 5 1 1 14 3.35770e+04
#> 6 1 2 14 2.94240e+04
#> 7 1 1 15 7.11000e+02
#> 8 1 2 15 4.80000e+01
#> 9 1 0 35 2.19321e-01
#> 10 1 0 36 2.87284e-01
The *_out.txt
file is usually the file you are
interested in, as it contains the measurement data for each survey you
defined in the input list. See
here for more information.
Read results
OpenMalaria’s output files use numerical indiccators in the second column which can make the interpretation of the values difficult for new users. Thus OMU provides a helper function which reads an output file and assigns column names and translates certain values (e.g. the survey dates).
outfile <- readOutputFile("exp_test/outputs/exp_test_base_out.txt")
head(outfile)
#> survey_date third_dimension measure value
#> 1: 1999-12-27 1 nHost 607
#> 2: 1999-12-27 2 nHost 3393
#> 3: 1999-12-27 1 nPatent 137
#> 4: 1999-12-27 2 nPatent 1155
#> 5: 1999-12-27 1 nUncomp 33577
#> 6: 1999-12-27 2 nUncomp 29424