search

PoonLab / twt

Generic framework for the coalescent simulation of pathogen trees within host trees in R
GNU Affero General Public License v3.0
3 stars 2 forks source link

sim.inner.tree() modifies its input #119

Open davidchampredon opened 4 years ago

davidchampredon commented 4 years ago

This is a follow-up of #118 where host "ID73" was the index case of the entire epidemic. My output from the ABM did not identify this host, as spotted by @ArtPoon , now the code is fixed.

I have noticed a strange behaviour, that is probably unintentional. It seems that sim.inner.tree() modifies its input. Below is the MRE. I would not expect the second call plot(outer.tree) to fail -- likely that outer.tree was changed on the way... Thanks!

library(twt)
set.seed(1234)
settings <- yaml.load_file('../data/sampled-hosts-1.yaml')
m     <- Model$new(settings)

outer.tree <- twt::load.outer.tree(model = m)
print('before')
plot(outer.tree)  # <-- OK
inner <- sim.inner.tree(outer.tree)
print('after')
plot(outer.tree)  # <-- fails
print('end')

and the associated YAML file:

InitialConditions:
  originTime: 3493
  size:
    host: 560
  indexType: 'host'

CompartmentTypes:
  'host':
    migration.rates: ()
    transition.rates: ()
    branching.rates: (host=1e-06)
    effective.size: 100
    generation.time: 1
    bottleneck.size: 1
    wait.time.distr:
      distr: 'exp'
      hyperparameters:
        - rate: 0.001
    popn.growth.dynamics:
      piece1:
        - startTime: 0
        - startPopn: 1
        - endTime: 5
        - endPopn: 20
      piece2:
        - startTime: 5
        - startPopn: 10

LineageTypes:
 'virus':

Compartments:
  'ID18':
    type: host
    source: 'ID73'
    branching.time: 2987
    unsampled: TRUE
  'ID277':
    type: host
    source: 'ID18'
    branching.time: 2982
    unsampled: TRUE
  'ID165':
    type: host
    source: 'ID18'
    branching.time: 2561
    unsampled: TRUE
  'ID221':
    type: host
    source: 'ID277'
    branching.time: 2556
    replicates: 1
  'ID190':
    type: host
    source: 'ID221'
    branching.time: 2129
    replicates: 1
  'ID105':
    type: host
    source: 'ID18'
    branching.time: 1564
    replicates: 1
  'ID188':
    type: host
    source: 'ID105'
    branching.time: 1553
    replicates: 1
  'ID406':
    type: host
    source: 'ID105'
    branching.time: 1549
    replicates: 1
  'ID192':
    type: host
    source: 'ID188'
    branching.time: 1542
    unsampled: TRUE
  'ID97':
    type: host
    source: 'ID188'
    branching.time: 1540
    unsampled: TRUE
  'ID207':
    type: host
    source: 'ID192'
    branching.time: 1534
    replicates: 1
  'ID399':
    type: host
    source: 'ID221'
    branching.time: 1532
    unsampled: TRUE
  'ID56':
    type: host
    source: 'ID192'
    branching.time: 1524
    replicates: 1
  'ID57':
    type: host
    source: 'ID56'
    branching.time: 1519
    unsampled: TRUE
  'ID218':
    type: host
    source: 'ID57'
    branching.time: 1511
    unsampled: TRUE
  'ID228':
    type: host
    source: 'ID218'
    branching.time: 1508
    replicates: 1
  'ID335':
    type: host
    source: 'ID218'
    branching.time: 1508
    unsampled: TRUE
  'ID149':
    type: host
    source: 'ID57'
    branching.time: 1505
    replicates: 1
  'ID235':
    type: host
    source: 'ID218'
    branching.time: 1504
    unsampled: TRUE
  'ID281':
    type: host
    source: 'ID235'
    branching.time: 1488
    unsampled: TRUE
  'ID288':
    type: host
    source: 'ID281'
    branching.time: 1478
    unsampled: TRUE
  'ID233':
    type: host
    source: 'ID288'
    branching.time: 1474
    unsampled: TRUE
  'ID152':
    type: host
    source: 'ID288'
    branching.time: 1459
    unsampled: TRUE
  'ID260':
    type: host
    source: 'ID288'
    branching.time: 1448
    unsampled: TRUE
  'ID418':
    type: host
    source: 'ID260'
    branching.time: 1436
    unsampled: TRUE
  'ID84':
    type: host
    source: 'ID105'
    branching.time: 1414
    unsampled: TRUE
  'ID237':
    type: host
    source: 'ID406'
    branching.time: 1319
    unsampled: TRUE
  'ID1':
    type: host
    source: 'ID233'
    branching.time: 1271
    replicates: 1
  'ID212':
    type: host
    source: 'ID277'
    branching.time: 1249
    unsampled: TRUE
  'ID463':
    type: host
    source: 'ID84'
    branching.time: 1229
    replicates: 1
  'ID389':
    type: host
    source: 'ID56'
    branching.time: 1227
    replicates: 1
  'ID93':
    type: host
    source: 'ID212'
    branching.time: 969
    unsampled: TRUE
  'ID372':
    type: host
    source: 'ID463'
    branching.time: 963
    unsampled: TRUE
  'ID214':
    type: host
    source: 'ID372'
    branching.time: 958
    unsampled: TRUE
  'ID474':
    type: host
    source: 'ID372'
    branching.time: 952
    unsampled: TRUE
  'ID295':
    type: host
    source: 'ID212'
    branching.time: 766
    unsampled: TRUE
  'ID100':
    type: host
    source: 'ID295'
    branching.time: 759
    unsampled: TRUE
  'ID21':
    type: host
    source: 'ID295'
    branching.time: 748
    replicates: 1
  'ID225':
    type: host
    source: 'ID295'
    branching.time: 742
    unsampled: TRUE
  'ID529':
    type: host
    source: 'ID474'
    branching.time: 593
    unsampled: TRUE
  'ID210':
    type: host
    source: 'ID73'
    branching.time: 592
    replicates: 1
  'ID426':
    type: host
    source: 'ID335'
    branching.time: 569
    replicates: 1
  'ID183':
    type: host
    source: 'ID277'
    branching.time: 528
    unsampled: TRUE
  'ID316':
    type: host
    source: 'ID183'
    branching.time: 519
    replicates: 1
  'ID61':
    type: host
    source: 'ID192'
    branching.time: 509
    unsampled: TRUE
  'ID553':
    type: host
    source: 'ID426'
    branching.time: 505
    replicates: 1
  'ID102':
    type: host
    source: 'ID57'
    branching.time: 415
    unsampled: TRUE
  'ID229':
    type: host
    source: 'ID102'
    branching.time: 410
    unsampled: TRUE
  'ID274':
    type: host
    source: 'ID229'
    branching.time: 401
    replicates: 1
  'ID191':
    type: host
    source: 'ID229'
    branching.time: 397
    unsampled: TRUE
  'ID104':
    type: host
    source: 'ID210'
    branching.time: 367
    unsampled: TRUE
  'ID430':
    type: host
    source: 'ID235'
    branching.time: 344
    replicates: 1
  'ID113':
    type: host
    source: 'ID57'
    branching.time: 337
    unsampled: TRUE
  'ID117':
    type: host
    source: 'ID430'
    branching.time: 325
    unsampled: TRUE
  'ID337':
    type: host
    source: 'ID152'
    branching.time: 279
    replicates: 1
  'ID33':
    type: host
    source: 'ID152'
    branching.time: 125
    unsampled: TRUE
  'ID560':
    type: host
    source: 'ID529'
    branching.time: 77
    unsampled: TRUE
  'ID73':
    type: host
    source: NA
    branching.time: NA
    unsampled: TRUE

Lineages:
  '1':
    type: virus
    sampling.time: 1263
    location: ID221
    replicates: 1
  '2':
    type: virus
    sampling.time: 483
    location: ID190
    replicates: 1
  '3':
    type: virus
    sampling.time: 1280
    location: ID105
    replicates: 1
  '4':
    type: virus
    sampling.time: 183
    location: ID188
    replicates: 1
  '5':
    type: virus
    sampling.time: 156
    location: ID406
    replicates: 1
  '6':
    type: virus
    sampling.time: 1467
    location: ID207
    replicates: 1
  '7':
    type: virus
    sampling.time: 1021
    location: ID56
    replicates: 1
  '8':
    type: virus
    sampling.time: 1486
    location: ID228
    replicates: 1
  '9':
    type: virus
    sampling.time: 1130
    location: ID149
    replicates: 1
  '10':
    type: virus
    sampling.time: 329
    location: ID1
    replicates: 1
  '11':
    type: virus
    sampling.time: 831
    location: ID463
    replicates: 1
  '12':
    type: virus
    sampling.time: 571
    location: ID389
    replicates: 1
  '13':
    type: virus
    sampling.time: 706
    location: ID21
    replicates: 1
  '14':
    type: virus
    sampling.time: 223
    location: ID210
    replicates: 1
  '15':
    type: virus
    sampling.time: 492
    location: ID426
    replicates: 1
  '16':
    type: virus
    sampling.time: 0
    location: ID316
    replicates: 1
  '17':
    type: virus
    sampling.time: 497
    location: ID553
    replicates: 1
  '18':
    type: virus
    sampling.time: 38
    location: ID274
    replicates: 1
  '19':
    type: virus
    sampling.time: 307
    location: ID430
    replicates: 1
  '20':
    type: virus
    sampling.time: 102
    location: ID337
    replicates: 1
ArtPoon commented 4 years ago

I couldn't reproduce this issue because a different exception gets thrown:

> settings <- yaml.load_file('~/git/twt/working/issue119.yaml')
> model <- Model$new(settings)
Error in FUN(X[[i]], ...) : 
  ID7316 of Lineage 16 is not a specified Compartment object.
Available compartments:
ID18_1ID277_1ID165_1ID221_1ID190_1ID105_1ID188_1ID406_1ID192_1ID97_1ID207_1ID399_1ID56_1ID57_1ID218_1ID228_1ID335_1ID149_1ID235_1ID281_1ID288_1ID233_1ID152_1ID260_1ID418_1ID84_1ID237_1ID1_1ID212_1ID463_1ID389_1ID93_1ID372_1ID214_1ID474_1ID295_1ID100_1ID21_1ID225_1ID529_1ID210_1ID426_1ID183_1ID316_1ID61_1ID553_1ID102_1ID229_1ID274_1ID191_1ID104_1ID430_1ID113_1ID117_1ID337_1ID33_1ID560_1ID73_1
ArtPoon commented 4 years ago

Never mind - this is a copy-paste error of the MRE

ArtPoon commented 4 years ago

Ran into a different problem:

> set.seed(1234)
> outer.tree <- load.outer.tree(model)
> plot(outer.tree)
Error in max(sapply(comp$get.lineages(), function(line) line$get.sampling.time())) : 
  invalid 'type' (list) of argument

This occurs because the comp object has no Lineages:

> node
[1] "ID316_1"
> comp
<Compartment>
  Public:
    add.lineage: function (new.lineage) 
    clone: function (deep = FALSE) 
    copy: function (deep = FALSE) 
    get.branching.time: function () 
    get.lineages: function () 
    get.name: function () 
    get.sampling.time: function () 
    get.source: function () 
    get.type: function () 
    initialize: function (name = NA, type = NA, source = NA, branching.time = NA, 
    is.unsampled: function () 
    remove.lineage: function (ex.lineage) 
    set.branching.time: function (new.branching.time) 
    set.sampling.time: function () 
    set.source: function (new.source) 
    set.type: function (new.type) 
    set.unsampled: function (is.unsampled) 
  Private:
    branching.time: 519
    deep_clone: function (name, value) 
    lineages: list
    name: ID316_1
    sampling.time: NA
    source: Compartment, R6
    type: CompartmentType, R6
    unsampled: FALSE
> comp$get.lineages()
list()
ArtPoon commented 4 years ago

Reproduced the problem:

> library(twt)
> set.seed(1234)
> #settings <- yaml.load_file('../data/sampled-hosts-1.yaml')
> settings <- yaml.load_file('issue119.yaml')
> m     <- Model$new(settings)
> 
> run <- twt::load.outer.tree(model = m)
> plot(run)  # <-- OK
> 
> events <- run$get.eventlog()$get.all.events()
> dim(events)
[1] 57  8
> 
> inner <- sim.inner.tree(run)
> 
> plot(run)  # <-- fails
Error in .plot.outer.tree(run, ...) : 
  Detected multiple roots in Run object: ID73_1, ID529_1
> 
> events2 <- run$get.eventlog()$get.events('transmission')
> dim(events2)
[1] 38  8

Very strange - sim.inner.tree has removed a bunch of the transmissions from the eventlog...

ArtPoon commented 4 years ago

@ewong347 ran into the same problem mentioned above in PoonLab/simclone#4 (private repo)

ewong347 commented 4 years ago

Another Yaml specification MRE


#Modelling from initial infection

InitialConditions:
  originTime: 10.0  # 
  size:
    Active: 1000  # number of susceptible at the origin time of the simulation for Type `active`
    latent: 100
  indexType: Active

CompartmentTypes:
  'Active':
    branching.rates: (Active=0.5, latent=0.05) 

    transition.rates: () #
    coalescent.rate: 1.5
    migration.rates: ()                
    generation.time: 1
    bottleneck.size: 1                 # bottleneck size when transmission event occurs
    effective.size : 2

  'latent':
    branching.rates: (Active=0.0005, latent=  0.00005)   
    transition.rates: ()
    migration.rates: ()                # no migration (assume no superinfection)
    generation.time: 1
    bottleneck.size: 1                # bottleneck size when transmission event occurs
    effective.size : 2

Compartments:
  'Active':
    type: Active                         # CompartmentType 
    replicates: 1
    source: NA
    branching.time: 0 
  'Latent':
    type: latent
    replicates: 1
    source: NA #try dropping me
    branching.time: 0 

Lineages:
  'A':
    type: virus
    sampling.time: 0          
    location: Active
    replicates: 1                   # Starts in the Active compartment
ArtPoon commented 4 years ago

I think this is related to #92. We should probably be cloning the Run object produced by sim.outer.tree before it is modified by sim.inner.tree which each call of the latter.