VLStm.jl observations

[t-pollington]

Hi @LloydChapman,

Please pardon the brevity:

CalcHHDists()

• Is abs() necessary here as it doesn't appear for the Haversine?

KnlHH()

• Is ones(nHH,nHH) correct here as shouldn't the diagonal contain zeros rather than ones, since the only way that intra-hhld transmission should be represented is in the \delta*1{ij}_ term. This will mean that KHH_{ij} will be non-zero for i=j when it should have been zero. Also K0 will be different. Actually this would apply to all typs and a second line of code should set the diagonals of KHH to zero.
• In K0 calc (see SI p3of39) shouldn't sumsumK=n be sumsumK=n^2 or n(n-1) instead? As otherwise an analysis of 3 paras vs 4 paras with the non-linear square relation of pairwise i,j comparisons vs linear increases in n may create unlike-for-like comparisons?

Tim.

[t-pollington]

Note also in newer Julia, the following (which appears twice) as CSV.read("data_final2.csv") should be replaced as CSV.read("data_final2.csv", DataFrame).

Also join() on a DataFrame no longer works, but assuming that innerjoin() is fine (when this occurs twice) in PlotSims.jl.

[t-pollington]

In functions.jl what does IM_IN do on line 102?

[t-pollington]

In functions.jl line 100 why is 8-status assigned to Event? ie why is Event treated like a counter in a relative sense when it could just be set to different values that describe what that event was?

[t-pollington]

In functions.jl line 103, why is only the first person sharing the internal migration of IM_IN used? Is this because multiple rows can account for the same person?

[t-pollington]

In tI[i] = t+rand(MersenneTwister(),NegativeBinomial(r1,p1))+1 why not allow infection to start at exposure? Is it better to use a truncation on the NB to restrict the support to >=1 instead, than doing this +1 correction?

[t-pollington]

"I think there's an issue here, as coded like this their PKDL onset could be after their death. Simulate birth and death times? Otherwise I'm assuming death rate is unaffected by KA." I agree. Also natural death from non-VL causes as a function of age and sex not included. Best for those who do get infected, for risk of death to also be modelled but only for the infectious part rather than the asyx part.

[t-pollington]

Why no +1 on this line tP[i] = tDor[i] + rand(MersenneTwister(),NegativeBinomial(r3,p3))

[t-pollington]

How did you obtain r0 and p0 since these are output as vectors from the following code:

    OTparams = readdlm(joinpath(@__DIR__,"OTparams.csv"),',',Float64)
    r0 = OTparams[:,1]
    p0 = OTparams[:,2]
    # r0 = 1.369639756715030*ones(9,1)
    # p0 = 0.380308479571719*ones(9,1)

Also 1.36 and 0.38 don't match with mean() or median() of the vector, nor the 1.34 value (0.38 was correct) in the SI just after Eqn S5.

[t-pollington]

🐛 In tI[i] = t + rand(MersenneTwister(),NegativeBinomial(r1,p1)) + 1 I think r1 and p1 are incorrect params for the exposure to onset time. Shouldn't this be r=3 and p=1.34ish instead?

[t-pollington]

tR[Asx] = t .+ rand(MersenneTwister(),Geometric(p2),A[1]) .+ 1
    tI[Exp] = t .+ randSizeBiasedNegativeBinomial(r1,p1,E[1])
    tR[Exp] = tI[Exp] + rand(MersenneTwister(),NegativeBinomial(r0[1],p0[1]),E[1]) .+ 1
    tDor[infDor] = t .+ randSizeBiasedNegativeBinomial(r0[1],p0[1],sum(infDor))
    tP[infDor] = tDor[infDor] + rand(MersenneTwister(),NegativeBinomial(r3,p3),sum(infDor))
    tR[infDor] = tP[infDor] + rand(MersenneTwister(),NegativeBinomial(r5,p5),sum(infDor)) .+ 1
    tR[infRec] = t .+ randSizeBiasedNegativeBinomial(r0[1],p0[1],sum(infRec))
    tP[Dor] = t .+ randSizeBiasedNegativeBinomial(r3,p3,D[1])
    tR[Dor] = tP[Dor] + rand(MersenneTwister(),NegativeBinomial(r5,p5),D[1]) .+ 1
    tR[PKDL] = t .+ randSizeBiasedNegativeBinomial(r5,p5,P[1])

Some missing +1. Is this OK?

[t-pollington]

🐛 Need to consider when at some times the tR[infRec] = t .+ randSizeBiasedNegativeBinomial(r0[1],p0[1],sum(infRec)) will evaluate to 0×1 Array{Int64,2} if randSizeBiasedNegativeBinomial(...,n=0). Better for these functions not to be performed and alternative control structures to guide the program away from using them, as 0×1 Array{Int64,2} may lead to instability when given as later function input.

[t-pollington]

🐛 As model can only work with integer months, randSizeBiasedNegativeBinomial() shouldn't accept non-integer r.

[t-pollington]

According to the parameters in the SI which are named as in Table S2 on p8, there is not a parameter for the success probability of the KA IP NB distribution. It is loaded in RunSims.jl as p1_str = "p1_DELTA_HGHR_PRESX_ASX_INFCTS_AllParas.csv" and is called p1 in the code, not to be confused with p1 in Table S2.

[t-pollington]

Why isn't r0[1] and p0[1] used when medians could be taken of these and used instead? It appears r0[1] & p0[1] are used initially for t=w=12 since r0[n] & p0[n] relate to the estimate in the study year n. Is this right? If so why is it that r0[9]=1 exactly and p0[9]=0.5 exactly, and why doesn't SI p5 section D mention the year-specific estimates for r0 & p0? Note these appear as r1 and p1 in the SI text.

[t-pollington]

Need to understand Event codes better and why they are used to advance individuals through Status using Status + Event.

[t-pollington]

Why is recovery from KA, and, AsxInf -> AsxDor grouped together for event 3?

elseif (inf[i] & infRec[i]) | (Asx[i] & AsxDor[i]) # add 3 for direct recovery from KA, or for progression from asx infection to dormant infection

since those who have recovered from KA are not infectious. Is it because AsxDor aren't infectious either but could become inf if they progress to PKDL?

[t-pollington]

I don't understand:

Event[i] = Status[IM_OUT[IM_IN.==i][1]] - Status[i]

[t-pollington]

elseif !(t<tIM[i]) & !(t>tD[i]) & !(t>tEM[i]) # individual didn't internally migrate in and hasn't already died/emigrated
            Event[i] = 1 # add 1 for birth, progression to KA, progression from KA to dormant infection, progression to PKDL, and recovery from PKDL
        end

Doesn't the above also add +1 if you are in Status=7 ie Rec and will move you to Status=8 ie death?