Antimalarial resistance/etf

[tbreweric]

Overview
I have been working to update malariasimulation to be able to simulate the effects of resistance to the antimalarial drugs used in clinical treatment. The updates to the model and code are integrating the modifications designed by Slater et al. (2016) (see https://malariajournal.biomedcentral.com/articles/10.1186/s12936-015-1075-7) without the need for additional infectious states/compartments. Briefly, the model for antimalarial resistance captures resistance to both the artemisinin and partner drug components of artemisinin combination therapies and five possible resistance outcomes:

1. Slow Parasite Clearance (artemisinin resistance)
2. Early Treatment Failure (artemisinin resistance)
3. Late Clinical Failure (partner drug resistance)
4. Late Parasitological Failure (partner drug resistance)
5. Reinfection During Prophylaxis (partner drug resistance)

The antimalarial_resistance/etf branch was set-up to develop the changes to malariasimulation required to capture early treatment failure (ETF). ETF causes individuals with artemisinin resistant infections to fail treatment and develop clinical disease. In malariasimulation terms, this represents people who have been selected to received drugs continuing into the D state rather than Tr (in the same way as already occurs as individuals fail treatment due to drug efficacy).

Rather than tracking each individuals infection status (non-resistant / resistant) in a new variable, we opted to apply this change in a simpler, probabilistic way in which each individual is assigned a probability of being successfully treated given the proportion of total infections that are resistant to artemisinin and the probability that artemisinin resistant infections result in early treatment failure. We can then use bernoulli trials (via bernoulli_multi_p()) to determine which individuals succeed/fail treatment. This additional behaviour was integrated into the existing calculate_treated() function (human_infection.R) prior to the drug efficacy success/failure steps.

In this PR I have:

1. Added a new script called antimalarial_resistance.R to house all new functions involved in the simulation of antimalarial resistance.
2. Written a function called set_antimalarial_resistance() (in antimalarial_resistance.R) that modifies the malariasimulation parameter list to append resistance parameters
3. Amended the calculate_treated() function (human_infection.R) to simulate the failure of treatment due to artemisinin resistance (early treatment failure, in which individuals destined for the Tr state fail treatment and instead contract clinical disease and go to D).
4. Written unit tests for set_antimalarial_resistance() in a new testthat script called test-antimalarial-resistance.R and for the updated calculate_treated() function (in test-infection-integration.R as this was where the existing test for the original function was saved).

antimalarial_resistance.R & set_antimalarial_resistance()
I have populated antimalarial_resistance.R with a single function, set_antimalarial_resistance(). set_antimalarial_resistance() is a user-facing function used to update the malariasimulation parameter list with resistance parameters for a given antimalarial drug (e.g. SP_AQ). Like set_clinical_treatment(), set_antimalarial_resistance() maps to drug parameters called by the set_drugs() function and is designed to be called once per drug that resistance is being parameterised for. The set_antimalarial_resistance() function accepts 10 inputs, as detailed in the Roxygen2 parameter definitions.

The function appends 7 resistance parameters. The artemisinin_resistance and partner_drug_resistance parameters describe the proportion of malaria infections in the given time step which are resistant to the artemisinin and partner drug component of the ACT (all drugs for which parameters are currently available in malariasimulation are ACTs) resistance is being parameterised for. The slow_parasite_clearance_prob and early_treatment_failure_prob parameters represent the probability that a case of artemisinin resistant malaria results in either slow parasite clearance or early treatment failure. The late_clinical_failure_prob, late_parasitological_prob, and reinfection_prob describe the probability that a case of partner drug resistant malaria results in the late clinical failure, late parasitological failure, or reinfection during prophylaxis respectively.

set_antimalarial_resistance() includes four checks which terminate the function if failed. The first checks that the number of values for each parameter input match the number of timesteps input (e.g. if timesteps specifies two updates to the resistance parameters, each parameter should have two values specifying the value to be updated to). The second and third check that the resistance proportions and probabilities all fall between 0 and 1. The final check ensures the drug index for which parameters are being assigned does not exceed the number of drugs parameterised for using set_drugs().

The function first turns antimalarial_resistance to TRUE. It then checks whether an antimalarial_resistance_drug index has already been assigned for the input drug and, if not, assigns it the next available integer. Once the index has been assigned, the input parameters are appended to the parameters list. The function concludes by outputting the amended parameters list which now has resistance parameters appended.

calculate_treated() updates
Early treatment failure is assumed to result in the immediate failure of treatment such that, in the malariasimulation model, individuals that receive clinical treatment fail this treatment and progress to the untreated clinical treatment infectious state (similar to failure due to drug efficacy).

The calculate_treated() function has been amended to enable individuals to fail clinical treatment due to resistance to the artemisinin resistant component of the drug they have received. The updated version has a new section that assigns resistance parameters from the parameter list given the drug that an individual has received. The artemisinin resistance and early treatment probability for the drug each individual has received are multiplied to get the probability that the treatment fails due to early treatment failure. This probability is subtracted from one to get the probability the infection is suceptible to treatment, and this is used to run bernoulli trials to determine who fails treatment. Those who remain susceptible then go through the existing drug-efficacy treatment success/failure steps. Those whose treatment is successful then have their infectious state, infectivity, drug, and drug time variables updated.

Unit Testing
I have created a new file in tests/testthat called test-antimalarial-resistance.R and added 6 tests. The first checks that the set_antimalarial_resistance() function turns the resistance toggle on in the parameter list. The second, third and fourth tests check that the in-built tests and by erroring out if the parameter inputs differ in length to the timesteps input, the resistance proportions fall outside the range 0-1, and or the resistance phenotype probabilities call outside the range 0-1. The fifth test checks that the function errors out if the input drug index exceeds the number of drugs for which parameters have been assigned using set_drugs(). The final test checks that the function assigns the parameters correctly when the order resistance parameters are assigned differs from the order they are assigned using set_drugs().

I have amended the test-infection-integration.R file to update the fourth test, ('calculate_treated correctly samples treated and updates the drug state') to account for the changes to the calculate_treated() function.

I have also added five additional tests for the updated calculate_treated() function. The first four check that the function returns an empty Bitset when there is no clinical treatment coverage, when the input Bitset of clinically_infected individuals is empty, when the input parameter list contains no resistance or clinical treatment parameters, or when all attempted treatment fails due to early treatment failure. The final test checks that the additional columns specified to be rendered by calculate_treated() are added to the output dataframe.

Outstanding Issues:
Some of the values rendered in the updated calculate_treatment() function (e.g. n_treat_eff_fail, n_treat_success) were designed for testing and will likely want to removed before integrating with dev.

The calculate_treated() function currently has resistance and treatment failure as separate calculations. However, the order of these will affect the number of individuals that fail treatment due to either, so we may want to combine the two to get a probability of failing treatment and run bernoulli_multi_p() with this combined probability and just output an n_treat_fail.

Note
This is the first time I've developed code for a collaborative project and submitted a PR so if you have recommendations on how I can improve / make the PR's easier for reviewers etc. please let me know!!

[giovannic]

Thank you for this PR! The comments and the testing is very thorough.

Could you set the PR destination to dev? That way github can run the tests and we can easily spot any glaring issues. We did make exceptions for the vivax and vignette work, but this seems more self contained.

[giovannic]

These lines are a bit long

[pwinskill]

Agree, this is tricky to parse. You call match_timestep(ts =parameters$antimalarial_resistance_timesteps[[drug_index[i]]], t = timestep) twice with the same parameters, so that can defo come out onto its own line.

I've a sneaking suspicion we can vectorise this whole lot (get rid of the for loop), but might be helpful to simplify first and then see if we can do that.

[tbreweric]

I've shortened these lines by pulling the match_timestep() step out which makes it more readable:

  matched_t <- match_timestep(ts = parameters$antimalarial_resistance_timesteps[[drug_index[i]]], t = timestep)
  artemisinin_resistance_proportion[drug_indices[[i]]] <- parameters$prop_artemisinin_resistant[[drug_index[i]]][matched_t]
  early_treatment_failure_probability[drug_indices[[i]]] <- parameters$early_treatment_failure_prob[[drug_index[i]]][matched_t]

I had another look at vectorising this step, which I agree should be possible, but couldn't puzzle it out yet.

[tbreweric]

I didn't merge to dev because there will be a series of additional PRs for antimalarial resistance for, at least, each of the resistance outcomes (slow parasite clearance, late clinical failure, late parasitological failure, and reinfection during prophylaxis. It seemed easier to me to merge them all into a single antimalarial_resistance branch, make sure they all talk to each other, and then merge that into dev when we're happy. However, if on reading this you still think it would be better to merge to dev I'll do that instead!

[giovannic]

we document outputs in the documentation for run_simulation

[tbreweric]

I've added documentation for the new outputs to the run_simulation() function in model.R. However, some/potentially all of these are likely to be removed as they were written for testing, so will need to remember to remove them once the antimalarial resistance work has been completed and we settle on a set of new outputs to include.

[tbreweric]

[giovannic]

I like how thorough this is, but it would be preferable to have deterministic unit tests. That way we can re-create and correct issues we may have in the future.

I would suggest picking one or two drug orderings which would convince us that set_antimalarial_resistance is not falling over. You can hard code different resistance parameters so that you can test that they are set in the correct order.

[tbreweric]

I've made the test deterministic and picked a single drug order/permutation:

test_that('set_antimalarial_resistance() assigns parameters correctly despite order in which resistance parameters are specified', {
  
  parameters <- get_parameters()
  parameters <- set_drugs(parameters = parameters, drugs = list(AL_params, SP_AQ_params, DHA_PQP_params))
  parameters <- set_clinical_treatment(parameters = parameters, drug = 2, timesteps = 1, coverages = 0.2)
  parameters <- set_clinical_treatment(parameters = parameters, drug = 1, timesteps = 1, coverages = 0.1)
  parameters <- set_clinical_treatment(parameters = parameters, drug = 3, timesteps = 1, coverages = 0.4)
  parameters <- set_antimalarial_resistance(parameters = parameters,
                                            drug = 2,
                                            timesteps = 1,
                                            artemisinin_resistance = 0.5,
                                            partner_drug_resistance = 0,
                                            slow_parasite_clearance_prob = 0.41,
                                            early_treatment_failure_prob = 0.2,
                                            late_clinical_failure_prob = 0,
                                            late_parasitological_prob = 0,
                                            reinfection_prob = 0)
  parameters <- set_antimalarial_resistance(parameters = parameters,
                                            drug = 3,
                                            timesteps = 1,
                                            artemisinin_resistance = 0,
                                            partner_drug_resistance = 0.43,
                                            slow_parasite_clearance_prob = 0,
                                            early_treatment_failure_prob = 0,
                                            late_clinical_failure_prob = 0.01,
                                            late_parasitological_prob = 0.42,
                                            reinfection_prob = 0.89)
  parameters <- set_antimalarial_resistance(parameters = parameters,
                                            drug = 1,
                                            timesteps = 1,
                                            artemisinin_resistance = 0.27,
                                            partner_drug_resistance = 0.61,
                                            slow_parasite_clearance_prob = 0.23,
                                            early_treatment_failure_prob = 0.9,
                                            late_clinical_failure_prob = 0.49,
                                            late_parasitological_prob = 0.81,
                                            reinfection_prob = 0.009)
  
  expect_identical(parameters$antimalarial_resistance, TRUE)
  expect_identical(unlist(parameters$antimalarial_resistance_drug), c(2, 3, 1))
  expect_identical(unlist(parameters$antimalarial_resistance_timesteps), rep(1, 3))
  expect_identical(unlist(parameters$prop_artemisinin_resistant), c(0.5, 0, 0.27))
  expect_identical(unlist(parameters$prop_partner_drug_resistant), c(0, 0.43, 0.61))
  expect_identical(unlist(parameters$slow_parasite_clearance_prob), c(0.41, 0, 0.23))
  expect_identical(unlist(parameters$early_treatment_failure_prob), c(0.2, 0, 0.9))
  expect_identical(unlist(parameters$late_clinical_failure_prob), c(0, 0.01, 0.49))
  expect_identical(unlist(parameters$late_parasitological_failure_prob), c(0, 0.42, 0.81))
  expect_identical(unlist(parameters$reinfection_during_prophylaxis), c(0, 0.89, 0.009))
  
})

I've also updated the test description to be more accurate of what is being tested.

[giovannic]

It is best to keep existing tests relatively unchanged so that it's clear to the reviewers that the new changes won't affect existing behaviour.

Could you revert this test, copy it into a new test, rename it and then test resistance behaviour in that new test.

[tbreweric]

I have reverted the test to it's original version (calculate_treated correctly samples treated and updates the drug state) and added a separate version of the same test for the antimalarial resistance behaviour ('calculate_treated correctly samples treated and updates the drug state when resistance set').

[giovannic]

It seems confusing to check both of these things when the first check is enough

Suggested change

	expect_identical(object = treated$to_vector(), expected = numeric(0), info = "Error: calculate_treated() returning non-zero number of treated individuals
	in the absence of clinical treatment")

[tbreweric]

I've removed the to_vector() test from test_that('calculate_treated returns empty Bitset when there is no clinical treatment coverage' ...) and kept the following:

expect_identical(object = treated$size(), expected = 0, info = "Error: calculate_treated() returning non-zero number of treated individuals in the absence of clinical treatment")

[giovannic]

Suggested change

	expect_identical(object = treated$to_vector(), expected = numeric(0), info = "Error: calculate_treated() returning non-zero number of treated individuals
	in the absence of clinically infected individuals")

[tbreweric]

I've removed the unnecessary additional test and left the following test:

expect_identical(object = treated$size(), expected = 0, info = "Error: calculate_treated() returning non-zero number of treated individuals in the absence of clinically infected individuals")

[giovannic]

Suggested change

	expect_identical(object = treated$to_vector(), expected = numeric(0), info = "Error: calculate_treated() returning non-zero number of treated individuals
	in the absence of clinical treatment or resistance parameters")

[tbreweric]

I've removed the unnecessary additional test and left the following test:

expect_identical(object = treated$size(), expected = 0, info = "Error: calculate_treated() returning non-zero number of treated individuals in the absence of clinical treatment or resistance parameters")

[giovannic]

Suggested change

	expect_identical(renderer$to_dataframe()[timestep,'n_ETF'], renderer$to_dataframe()[timestep,'n_treated'], info = "Error: Number of
	early treatment failures does not match number of treated individuals when artemisinin resistance proportion and
	and early treatment failure probability both equal 1")

[tbreweric]

Updated the test desc to "Number of treatment failures matches number of individuals treated when artemisinin resistance proportion and early treatment failure probability both set to 1" and implemented the suggested change to the test.

[giovannic]

It would be nice to have a vignette, but it could be done at a later date. We shouldn't add empty vignettes to the repo

[pwinskill]

+1, lets get this in with the PR, otherise as it isn't the "fun" bit it'll get forgotten!
It's also helpful for the reviewer to see a worked example in a vignette

[tbreweric]

I've written a draft vignette for the antimalarial resistance feature that covers how to use the set_antimalarial_resistance() function to parameterise runs with resistance and gives examples using the early treatment failure functionality currently built. This may need updating as additional resistance outcomes (e.g. a partner drug resistant outcome like late clinical failure) are added, but should give us something to go with!

[pwinskill]

This is looking great 🌈
A few addtional comments/suggestions/requests from me.
Thanks Tom!

[pwinskill]

I we have a single timesteps argument for both resistance proportion and pheotype prob, I think this always has to be an "and" right? Ie, as your check on line 26 shows, both must be equal in length to timesteps.

[tbreweric]

Yes you're right - I've changed this to "and"

[pwinskill]

Suggested change

	stop("Number of resistance parameter vectors do not match time steps specified for update")
	stop("Length of one or more resistance parameter vectors does not match time steps specified for update")

[tbreweric]

This reads better, I've made the change as recommended

[pwinskill]

Suggested change

	# Check the drug_index for the we're drug setting parameters for:
	# Check the drug_index for the drug we're setting parameters for:

[tbreweric]

This also reads better and I've also made the change as recommended

[pwinskill]

Is there a reason why parameters$antimalarial_resistance isn't in the deafault parameter list and set to FALSE?

[tbreweric]

I've added the antimalarial resistance parameters to the get_parameters() function with the following documentation and defaults:

[pwinskill]

Agree, this is tricky to parse. You call match_timestep(ts =parameters$antimalarial_resistance_timesteps[[drug_index[i]]], t = timestep) twice with the same parameters, so that can defo come out onto its own line.

I've a sneaking suspicion we can vectorise this whole lot (get rid of the for loop), but might be helpful to simplify first and then see if we can do that.

[pwinskill]

Suggested change

	susceptible_to_treatment_index <- bernoulli_multi_p(p = 1 - (artemisinin_resistance_proportion * early_treatment_failure_probability))
	unsuccessful_treatment_probability <- artemisinin_resistance_proportion * early_treatment_failure_probability
	susceptible_to_treatment_index <- bernoulli_multi_p(p = 1 - nsuccessful_treatment_probability)

[pwinskill]

Might be clearer?

[tbreweric]

Agreed, I've simplified as recommended

[pwinskill]

ETF will be unknown to most people. I'd probably call this n_early_treatment_failure to be very clear

[tbreweric]

Agreed, n_ETF changed to n_early_treatment_failure both as an object in calculate_treated() and in the run_simulation() outputs

[pwinskill]

Pull this out into a separate line with a helpful name?

[tbreweric]

Stored the drugs taken by those successfully treated in the following way:

successfully_treated_drugs <- drugs[successfully_treated_index]

and updated the queue_update() steps as follows:

variables$infectivity$queue_update(
  parameters$cd * parameters$drug_rel_c[successfully_treated_drugs],
  successfully_treated
)
variables$drug$queue_update(
  successfully_treated_drugs,
  successfully_treated
)

[pwinskill]

In general, I'd explicitly test for the error message you're expecting (like here)

Otherwise its easy to past the test as a result of a differnt (and potentially unexpected/wrong) error

[tbreweric]

I've changed the test to so that the error message matches that which is being tested in the set_antimalarial_resistance("Artemisinin and partner-drug resistance proportions must fall between 0 and 1") function:

"Artemisinin and partner-drug resistance proportions must fall between 0 and 1"

[pwinskill]

+1, lets get this in with the PR, otherise as it isn't the "fun" bit it'll get forgotten!
It's also helpful for the reviewer to see a worked example in a vignette

[tbreweric]

I've responded to each of the comments and implemented most/all of the updated changes now. Would you be able to review the changes and let me know whether I can merge this PR?

[giovannic]

Thank you so much for these changes! And the tests are quite extensive and clear, it's a joy to read.

Sorry to drag this on, I just want to be sure with calculate_treated, comments inline, and a few picky bits on documentation/outputs. Feel free to push back on the latter and correct me on the former if I'm mistaken.

[giovannic]

Suggested change

	if(parameters$antimalarial_resistance == TRUE) {
	if(parameters$antimalarial_resistance) {

[tbreweric]

Thanks for this additional set of comments @giovannic!

I have implemented the change here as recommended.

[giovannic]

I think there is a subtle bug here.

Is drug_indices supposed to be the index of the population who have been treated with drug parameters$antimalarial_resistance_drug[[i]] or i? Because I think it should be the former.

You can test this by copying your test called 'calculate_treated correctly samples treated and updates the drug state when resistance set', and adapting it to only set resistance on the drug at index 2. Hopefully that should reveal the bug and test any fix you make.

I've suggested a fix below and renamed some variables to clarify my thought process. I don't know if it would work without seeing the test result.

[Note] I think it's better to assign drug_indices and drug_index inside the loop since they are never used outside of the loop.

[Note 2] We could save on a which call for drug_indices and just use it as a boolean index with the same outcome

[Note 3] Am I right in thinking drug_index is just i?

Suggested change

	drug_indices <- list(); drug_index <- vector()
	for(i in seq_along(parameters$antimalarial_resistance_drug)) {
	drug_indices[[i]] <- which(drugs == i)
	drug_index[i] <- which(parameters$antimalarial_resistance_drug == i)
	matched_t <- match_timestep(ts = parameters$antimalarial_resistance_timesteps[[drug_index[i]]], t = timestep)
	artemisinin_resistance_proportion[drug_indices[[i]]] <- parameters$prop_artemisinin_resistant[[drug_index[i]]][matched_t]
	early_treatment_failure_probability[drug_indices[[i]]] <- parameters$early_treatment_failure_prob[[drug_index[i]]][matched_t]
	}
	for(i in seq_along(parameters$antimalarial_resistance_drug)) {
	drug <- parameters$antimalarial_resistance_drug[[i]]
	treated_with_drug <- drugs == drug
	matched_t <- match_timestep(
	ts = parameters$antimalarial_resistance_timesteps[[i]],
	t = timestep
	)
	artemisinin_resistance_proportion[treated_with_drug] <- parameters$prop_artemisinin_resistant[[i]][matched_t]
	early_treatment_failure_probability[treated_with_drug] <- parameters$early_treatment_failure_prob[[i]][matched_t]
	}

[tbreweric]

Great spot!

The purpose of this line of code is to, for each drug in drugs, retrieve the antimalarial resistance parameters associated with that drug in that timestep.

I have replaced the for-loop that retrieves the resistance parameters for each drug assigned to individuals receiving treatment with a function called get_antimalarial_resistance_parameters() (antimalarial_resistance.R). The function takes as it's inputs the parameters list, the vector of drugs, and the current timestep, and returns a list of vectors, each of length length(drugs), of each of the antimalarial resistance parameters corresponding to the drugs administered and the current timestep. To avoid issues arising when users do not assign resistance to all drugs for which clinical treatment is specified, within get_antimalarial_resistance_parameters() the vectors are initiated with no-resistance defaults (0's for the resistance proportions / probabilities, parameters$dt for the slow parasite clearance dt). The function then loops through each of the drugs for which resistance has been parameterised (parameters$antimalarial_resistance_drugs), gets the indices from drugs which matches the drug, and uses these to assign the appropriate resistance parameters to the output vectors.

To support the use of get_antimalarial_resistance_parameters(), which was pulled from the antimalarial_resistance/spc branch, I have also updated the set_antimalarial_resistance() function to include slow_parasite_clearance_time as an input.

I have added some additional tests for the get_antimalarial_resistance_parameters() function in test-antimalarial-resistance.R and updated those that existed for set_antimalarial_resistance() so that they are compatible with the updated version of the function. I have also gone through the unit tests I'd written in test-infection-integration.R to make them compatible with these changes. I have also added the additional test recommended above, in which resistance is only set to one of multiple drugs to test-infection-integration.R. Finally, I have updated the documentation in parameters.R for get_parameters() to include an entry for dt_slow_parasite_clearance and have updated the vignette (Antimalarial_Resistance.Rmd) to explain the new dt_slow_parasite_clearance parameter and ensure it's compatible with the new version of set_antimalarial_resistance()

Here is an example of some simulations run in which, for all 4, three ACTs are parameterised for clinical_treatment. The "none" line represents the scenario in which resistance isn't parameterised for any of the drugs, the "Single" when resistance is only parameterised for one, etc. For each drug, resistance is introduced (when parameterised) or the year corresponding to their drug index (e.g. drug 2 gets resistance in year two) - black geom_vlines() show the points where resistance is introduced. I think this demonstrates that the changes made avoid the bug where resistance parameter retrieval fails when the user does not call set_antimalarial_resistance() for all drugs for which clinical treatment has been specified.

[giovannic]

Fantastic!

[giovannic]

To get the automated checks to pass!

Suggested change

	# Load the requisite packages:
	#library(malariasimulation)
	# Use devtools to load the updated malariasimulation:
	devtools::load_all("C:/Users/trb216/OneDrive - Imperial College London/Desktop/malariasimulation/malariasimulation.Rproj")
	# Load the requisite packages:
	library(malariasimulation)

[tbreweric]

Good spot - removed the load_all() call, uncommented library(malariasimulation) and confirmed vignette can build!

[pwinskill]

Hey Tom,
Thanks for this - both code and PR were very clear for me to follow!
A few last suggestions/qs - nothing major.
🎈🍕

[pwinskill]

Is it ok to be adding this to dev, inferring that the user can set things like reinfection prob, when they actually have no impact?
Not sure best practice here, should they be removed or have warnings?

[tbreweric]

My plan was to merge to feat/antimalarial_resistance, develop each of the five different resistance outcomes individually on their own branches from feat/antimalarial_resistance, consolidate them on feat_antimalarial_resistance, and then merge to dev.

I initially opted to do this to because I was new to the model & GitHub, but now that I'm more comfortable with both I'm happy to do something more efficient . So if it is better to merge this branch to dev, and then branch from/merge back to dev to build in each additional resistance outcome then I'm happy to do that instead!

[giovannic]

I suggested we merge each outcome into dev directly. Since this seems to work in isolation, we can put it through automated tests earlier, we will have smaller PRs (with dev) and won't fall behind dev waiting for all 5 to be implemented.

And I agree, we should have an error if someone tries to use an unimplemented feature!

[pwinskill]

Suggested change

	artemisinin_resistance_proportion <- numeric(length = length(drugs))
	partner_drug_resistance_proportion <- numeric(length = length(drugs))
	slow_parasite_clearance_probability <- numeric(length = length(drugs))
	early_treatment_failure_probability <- numeric(length = length(drugs))
	late_clinical_failure_probability <- numeric(length = length(drugs))
	late_parasitological_failure_probability <- numeric(length = length(drugs))
	reinfection_during_prophylaxis_probability <- numeric(length = length(drugs))
	dt_slow_parasite_clearance <- rep(parameters$dt, length = length(drugs))
	empty_vector <- numeric(length = length(drugs))
	artemisinin_resistance_proportion <- empty_vector
	partner_drug_resistance_proportion <- empty_vector
	slow_parasite_clearance_probability <- empty_vector
	early_treatment_failure_probability <- empty_vector
	late_clinical_failure_probability <- empty_vector
	late_parasitological_failure_probability <- empty_vector
	reinfection_during_prophylaxis_probability <- empty_vector
	dt_slow_parasite_clearance <- rep(parameters$dt, length = length(drugs))

This is more effcient. Maybe worth changing if this function will be called often

[tbreweric]

This function will be called once per timestep for each timestep clinical treatment is used and resistance is parameterised. I quickly compared the original and suggested methods and it was a good bit faster so have implemented as recommended!

[pwinskill]

This check seems like overkill to me?
You specifiy that those vectors are length(drugs) by definition when you create them. Is there reason why they should/would change?

[tbreweric]

I don't think they should/would, but I added the test as insurance. You're probably right though, so have removed this check from the function.

[pwinskill]

Suggested change

	renderer$render('n_early_treatment_failure', n_early_treatment_failure, timestep)
	} else {
	susceptible_to_treatment <- seek_treatment
	n_early_treatment_failure <- n_treat - susceptible_to_treatment$size()
	renderer$render('n_early_treatment_failure', n_early_treatment_failure, timestep)
	}
	} else {
	susceptible_to_treatment <- seek_treatment
	n_early_treatment_failure <- n_treat - susceptible_to_treatment$size()
	}
	renderer$render('n_early_treatment_failure', n_early_treatment_failure, timestep)

Can you safely pull this render outside the ifelse here to avoid a duplicated line?

[tbreweric]

Yes - have pulled the lines calculating and rendering n_early_treatment_failure outside of the ifelse statement.

[pwinskill]

Looks good Tom 🦄
You'll need to tidy the merge conflict in the antimalarial_reisistance.R to pass checks

[tbreweric]

Finally wrangled all of the coding cats so that the checks have passed - are well still planning on merging to dev? Or do we want to merge to antimalarial_resistance so I can combine the ETF and SPC work first?