31:T505,// [Note] There is support in this notebook for setting a certain portion of the population as workers and assigning them workplaces, but that isn't currently being used since the current version of the game includes going to school and work as a subset of "large gatherings." ageDistribution = [{range: '0-19', portionOfPopInThisRange: 0.26, susceptibility: 0.003, asymptomaticPortion: 0.50 , hospitalIfSymptomatic: 0.02, employedProportion: 0.95}, // employed proportion includes kids going to 'school' -- though it's not currently being used in the game. {range: '20-39', portionOfPopInThisRange: 0.27, susceptibility: 0.008, asymptomaticPortion: 0.4 , hospitalIfSymptomatic: 0.16, employedProportion: 0.8}, {range: '40-59', portionOfPopInThisRange: 0.26, susceptibility: 0.010, asymptomaticPortion: 0.3 , hospitalIfSymptomatic: 0.25, employedProportion: 0.8}, {range: '60-79', portionOfPopInThisRange: 0.17, susceptibility: 0.025, asymptomaticPortion: 0.2, hospitalIfSymptomatic: 0.37, employedProportion: 0.5}, {range: '80+', portionOfPopInThisRange: 0.04, susceptibility: 0.035, asymptomaticPortion: 0.15, hospitalIfSymptomatic: 0.50, employedProportion: 0.05} ] 32:T45c,infectors = { let infectors= [] for (let ageGroup = 0 ; ageGroup < ageDistribution.length ; ageGroup++ ){ let susceptibilityWeightedProportion = ageDistribution[ageGroup].portionOfPopInThisRange * ageDistribution[ageGroup].susceptibility / d3.sum(ageDistribution, d => d.portionOfPopInThisRange * d.susceptibility) let peopleNeeded = roundProb(numInfectors * susceptibilityWeightedProportion) let asymptomaticNeeded = Math.round(peopleNeeded * ageDistribution[ageGroup].asymptomaticPortion) while (peopleNeeded > 0) { infectors.push({age: ageDistribution[ageGroup].range, asymptomatic: (asymptomaticNeeded > 0) ? 1 : 0 , hospitalized: (Math.random() { for (let ageGroup = 0 ; ageGroup < ageDistribution.length ; ageGroup++ ){ let peopleNeeded = roundProb(d.n * ageDistribution[ageGroup].portionOfPopInThisRange) while (peopleNeeded > 0) { susceptibles.push({age: ageDistribution[ageGroup].range, susceptibility: ageDistribution[ageGroup].susceptibility, worker: Math.random() < ageDistribution[ageGroup].employedProportion, housemate: 0, // this updates below - after inNeighborhood is done populating coworker: 0, // this updates below - after inNeighborhood is done populating neighborhood: d.type } ) peopleNeeded--; } } // assign housemates and coworkers from within the neighborhood if(d.type == 'inNeighborhood'){ _.sampleSize(_.keys(susceptibles), infector.nPeopleInHouse - 1).forEach(i => { susceptibles[i].housemate = 1 }) _.sampleSize(_.keys(susceptibles), roundProb(susceptibles.length / infectorDestinationsInNeighborhood)).forEach(i => { if(susceptibles[i].worker){susceptibles[i].coworker = 1} }) } }) let housemates = _.keys(_.pickBy(susceptibles, {housemate: 1})) let coworkers = _.keys(_.pickBy(susceptibles, {coworker: 1})) let inNeighborhood = _.keys(_.pickBy(susceptibles, {neighborhood: "inNeighborhood"})) let visitingInfector = _.keys(_.pickBy(susceptibles, {neighborhood: "visitingInfector"})) let infectorVisiting = _.keys(_.pickBy(susceptibles, {neighborhood: "infectorVisiting"})) // create object where all exposures will be stored. let exposures = {} // go through days of infectiousness and add to exposures Object. for(let dayNum = 1; dayNum <= daysWhileInfectious; dayNum++){ if(consoleLogging){console.log("Day: " + dayNum)} // set infectiousness to 1 for symptomatic, and to asymptomaticInfectiousness for a/pre-symptomatic let infectiousness = ( dayNum <= daysPresymptomaticInfectious || infector.asymptomatic ) ? asymptomaticInfectiousness : 1; // infect people in house if(infector.nPeopleInHouse > 1){ housemates.forEach(i => exposures[+i] = (exposures[+i] || 0) + (infectiousness * dailyHouseholdExposures)) } // now iterate through all 'go-out' decision types. Object.keys(decisionsPerWeek).forEach(type => { let baselineOutProb = decisionsPerWeek[type] / 7 if (type == 'work'){ // special case "work" if (infector.worker == 1 && Math.random() < baselineOutProb) { // if worker goes to work that day // grab proportion of coworkers also working that day and add to their exposures. //console.log(_.sampleSize(coworkers, roundProb(coworkers.length * baselineOutProb))) _.sampleSize(coworkers, roundProb(coworkers.length * baselineOutProb)).forEach(i => exposures[+i] = (exposures[+i] || 0) + infectiousness) // TO DO : need to also infect some number of other people who go to that location too. } else { if(consoleLogging){console.log("Not going out for work on day: " + dayNum)}} } else { let goOut = ( Math.random() < baselineOutProb && (infector.asymptomatic || dayNum <= daysPresymptomaticInfectious || (infector.hospitalized && dayNum <= (daysPresymptomaticInfectious + daysSymtpomaticInfectious - daysHospitalized) && Math.random() <= ChanceOfGoingOutIfSymptomatic ) || (infector.hospitalized === 0 && Math.random() <= ChanceOfGoingOutIfSymptomatic ) ) ) if(goOut){ if(consoleLogging){console.log(type + ': going out')} [...Array(type == 'large gathering' ? numLargeGatheringDestinations() : 1)].forEach(() => { let inNeighborhood = Math.random() < (infectorDestinationsInNeighborhood / (infectorDestinationsInNeighborhood + infectorDestinationsInOtherNeighborhoods)) // destinationMultipliers if(inNeighborhood){ _.concat( _.sampleSize(inNeighborhood, roundProb(inNeighborhood.length * destinationMultipliers.homeNeighborhood)), _.sampleSize(visitingInfector, roundProb(visitingInfector.length * destinationMultipliers.visitingNeighborhood)) ) .forEach(i => exposures[i] = (exposures[i] || 0) + (exposureByType[type] * infectiousness)) } else { _.concat( _.sampleSize(infectorVisiting, roundProb(infectorVisiting.length * destinationMultipliers.homeNeighborhood)), _.sampleSize(inNeighborhood, roundProb(inNeighborhood.length * destinationMultipliers.visitingNeighborhood)) ) .forEach(i => exposures[i] = (exposures[i] || 0) + (exposureByType[type] * infectiousness)) } if(consoleLogging){console.log(exposures)} }) } } }) } let cumulativeInfections = 0 let householdCumulativeInfections = 0 let coworkerCumulativeInfections = 0 Object.keys(exposures).forEach( d => { cumulativeInfections += (1 - Math.pow( ( 1 - susceptibles[d].susceptibility ), exposures[d]) ) if(susceptibles[d].housemate){ householdCumulativeInfections += (1 - Math.pow( ( 1 - susceptibles[d].susceptibility ), exposures[d]) ) } if(susceptibles[d].coworker){ coworkerCumulativeInfections += (1 - Math.pow( ( 1 - susceptibles[d].susceptibility ), exposures[d]) ) } }) //console.log(exposures) if(consoleLogging){ console.log("Infector:") console.log(infector) console.log("Susceptibles:") console.log(susceptibles) console.log('END OF INFECTOR. Infected ' + cumulativeInfections) } return { cumulativeInfections: cumulativeInfections, householdCumulativeInfections: householdCumulativeInfections, // not currently using separate workplaces in actual game so commenting this out // coworkerCumulativeInfections: coworkerCumulativeInfections, // nCoworkers: coworkers.length } } 34:[] 2f:{"title":"R0 Estimate for Infectious Disease Simulation","cells":[{"id":0,"value":"md`# R0 Estimate for Infectious Disease Simulation\nThis R0 Calculator was built to help tune the simulation used in the [*People of the Pandemic*](https://peopleofthepandemicgame.com/) game. (In depth description of methodology is [here](https://observablehq.com/@osserman/people-of-the-pandemic-simulation-methodology))\n\nSince the game doesn't track each subsequent infection from each infectious individual (mainly for computational / performance reasons), there is no direct data from the game that makes it easy to calcuate the R0 of the game -- not to mention other parameters like the Secondary Household Attack Rate.\n\nThis notebook takes a different approach, looking at 1K independent infectious individuals, each of whom has a completely susceptible population in which they operate, using the same logic (contact structure, social mixing structure, householding, disease course and parameters, etc) that the game uses in its no-mitigation / baseline scenario.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":274,"value":"html `

Estimate for the game's R0 using parameters below is ${R0.toFixed(2)}

\nThis is R0 for \"business as usual\" in the simulation assuming no mitigation.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":283,"value":"md`## Change parameters below to see how it impacts the R0`","pinned":false,"hidden":false,"mode":"ojs"},{"id":2,"value":"$31","pinned":true,"hidden":false,"mode":"ojs"},{"id":78,"value":"md`### Global Destination Numbers`","pinned":false,"hidden":false,"mode":"ojs"},{"id":51,"value":"viewof destinationsInNeighborhood= slider({\n min: 0, \n max: 15, \n step: 1, \n value: 7, \n title: 'Number of destination in each \"neighborhood\" within the zipcode.'\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":57,"value":"viewof peoplePerDestination= slider({\n min: 0, \n max: 100, \n step: 5, \n value: 20, \n title: 'Number of people in the neighborhood for each destination.'\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":53,"value":"viewof otherNeighborhoodsVisited= slider({\n min: 0, \n max: 10, \n step: 1, \n value: 3, \n title: 'How many neighborhoods beyond their own might people go to in the simulation.'\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":100,"value":"susceptibleInOtherNeighborhoods = otherNeighborhoodsVisited * destinationsInNeighborhood * peoplePerDestination","pinned":true,"hidden":false,"mode":"ojs"},{"id":403,"value":"susceptibleVisitingInfectorsNeighborhood = susceptibleInOtherNeighborhoods","pinned":true,"hidden":false,"mode":"ojs"},{"id":61,"value":"susceptibleInNeighborhood = destinationsInNeighborhood * peoplePerDestination - 1 // -1 to remove the infector from the count","pinned":true,"hidden":false,"mode":"ojs"},{"id":81,"value":"md`### Potential Destinations for individual infectors`","pinned":false,"hidden":false,"mode":"ojs"},{"id":72,"value":"viewof infectorDestinationsInNeighborhood= slider({\n min: 0, \n max: destinationsInNeighborhood, \n step: 1, \n value: 7, \n title: 'How many of the ' + destinationsInNeighborhood + ' destinations in their own neighborhood do they visit. '\n}) ","pinned":false,"hidden":false,"mode":"ojs"},{"id":75,"value":"viewof infectorDestinationsInOtherNeighborhoods= slider({\n min: 0, \n max: destinationsInNeighborhood * otherNeighborhoodsVisited, \n step: 1, \n value: 5, \n title: 'How many of the ' + destinationsInNeighborhood * otherNeighborhoodsVisited + ' destinations in the other ' + otherNeighborhoodsVisited + ' neighborhood do they visit. '\n}) ","pinned":false,"hidden":false,"mode":"ojs"},{"id":481,"value":"md `The next two fields describe the different activities people might do. First how many times they do each per week, and then the exposure intensity associated with an exposure from each of these activities.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":402,"value":"decisionsPerWeek = ({'food':3,'exercise':5, 'small gathering':2, 'large gathering':5, 'work':0})","pinned":true,"hidden":false,"mode":"ojs"},{"id":401,"value":"exposureByType = ({'food':0.3, 'exercise':0.1, 'small gathering':1, 'large gathering':0.3, 'work':1.5})","pinned":true,"hidden":false,"mode":"ojs"},{"id":118,"value":"viewof ChanceOfGoingOutIfSymptomatic= slider({\n min: 0, \n max: 1, \n step: 0.05, \n value: 0.5, \n title: \"Chance of going out if Symptomatic, but not in need of hospital care.\"\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":400,"value":"md` This formula is used to pick \"destinations\" for large gatherings: we model the decision to go to for large gatherings as going to 1-5 regular destinations in one day. `","pinned":false,"hidden":false,"mode":"ojs"},{"id":399,"value":"numLargeGatheringDestinations = function() {return _.random(1, 5)}","pinned":true,"hidden":false,"mode":"ojs"},{"id":332,"value":"md` This is formula used to assign number of people to houses. Currently it's just a random number between 1 and 5. `","pinned":false,"hidden":false,"mode":"ojs"},{"id":38,"value":"numPeopleInHouseFormula = function(){ return _.random(1, 5)}","pinned":true,"hidden":false,"mode":"ojs"},{"id":291,"value":"md`### Disease Parameters`","pinned":true,"hidden":false,"mode":"ojs"},{"id":18,"value":"viewof daysPresymptomaticInfectious = slider({\n min: 1, \n max: 5, \n step: 1, \n value: 2, \n title: \"Days of Presymptomatic Infectiousness.\",\n description: \"How many days someone is infectious on average before they develop symptoms\"\n})\n","pinned":false,"hidden":false,"mode":"ojs"},{"id":21,"value":"viewof daysSymtpomaticInfectious = slider({\n min: 5, \n max: 14, \n step: 1, \n value: 9, \n title: \"Days of Symptomatic Infectiousness.\",\n description: \"How many days someone is infectious and symptomatic\"\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":238,"value":"viewof daysHospitalized = slider({\n min: 1, \n max: daysSymtpomaticInfectious, \n step: 1, \n value: daysSymtpomaticInfectious-1, \n title: \"How Many Days of Hospitalization\",\n description: \"When Hospitalization Is Required How Many Days Does it Last (is bound by number of days of symptomatic infectiosness)\"\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":397,"value":"viewof asymptomaticInfectiousness = slider({\n min: 0, \n max: 1, \n step: 0.05, \n value: 0.4, \n title: \"Infectiousness of A/Pre-symptomatic Individuals\",\n description: \"When someone becomes infectious but are not symptomatic how infectious are they as a proportion of the infectious of a symptomatic individual.\"\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":121,"value":"daysWhileInfectious = daysPresymptomaticInfectious+daysSymtpomaticInfectious","pinned":true,"hidden":false,"mode":"ojs"},{"id":155,"value":"viewof dailyHouseholdExposures = slider({\n min: 1, \n max: 10, \n step: 1, \n value: 2, \n title: \"Daily Household Exposures.\",\n description: \"Number of times an infected person 'exposes' their housemates per day. \"\n})\n","pinned":false,"hidden":false,"mode":"ojs"},{"id":311,"value":"md`## Calculations Are Below`","pinned":false,"hidden":false,"mode":"ojs"},{"id":10,"value":"md`We'll start out with 1000 independent infectors that are randomly sampled based on age distribution proportions and the susceptibility of different age-groups. This mainly impacts how consistent the numbers are from run to run, and 1k seems like a good balance between being quick to run and stable from run to run.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":47,"value":"numInfectors = 1000","pinned":true,"hidden":false,"mode":"ojs"},{"id":314,"value":"md` This is our array of infectors.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":9,"value":"$32","pinned":false,"hidden":false,"mode":"ojs"},{"id":393,"value":"md` This is the proportion of in-neighborhood people and visitors who will be at any given destination on any given day. It's used in the infectionsByInfector() function.` ","pinned":false,"hidden":false,"mode":"ojs"},{"id":392,"value":"destinationMultipliers = {\n let dailyNonWorkDestinations = (d3.sum(Object.values(decisionsPerWeek)) - decisionsPerWeek.work) / 7\n let inNeighborhoodDestinationChance = (infectorDestinationsInNeighborhood / (infectorDestinationsInNeighborhood + infectorDestinationsInOtherNeighborhoods))\n return ({\n homeNeighborhood: dailyNonWorkDestinations * (inNeighborhoodDestinationChance) / destinationsInNeighborhood,\n visitingNeighborhood: dailyNonWorkDestinations * (1-inNeighborhoodDestinationChance) / ( destinationsInNeighborhood * otherNeighborhoodsVisited)\n} )\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":318,"value":"md` This is the function to calculate expected infections from each infector -- assuming no other infectors are also exposing the same people.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":34,"value":"$33","pinned":false,"hidden":false,"mode":"ojs"},{"id":389,"value":"md`This is the stats array from running the above infectionsByInfector() function on all infectors.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":388,"value":"infectionStats = infectors.map(d => infectionsByInfector(d))","pinned":false,"hidden":false,"mode":"ojs"},{"id":321,"value":"md` This is array with expected infections for each independ infector.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":111,"value":"infectionsByEachInfector = infectionStats.map(d => d.cumulativeInfections)","pinned":false,"hidden":false,"mode":"ojs"},{"id":325,"value":"md` R0 is just the average accross all those infectors`","pinned":false,"hidden":false,"mode":"ojs"},{"id":178,"value":"R0 = d3.mean(infectionsByEachInfector)","pinned":false,"hidden":false,"mode":"ojs"},{"id":385,"value":"md `Secondary Attack Rates for Households = cumulative household infections / cumulative number of housemates.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":384,"value":"householdAttackRate = d3.sum(infectionStats, d => d.householdCumulativeInfections) / d3.sum(infectors, d => d.nPeopleInHouse-1) ","pinned":false,"hidden":false,"mode":"ojs"},{"id":382,"value":"md `R0 by infector age is cumulative infections for each infector age group / total number of infectors of that age. `","pinned":false,"hidden":false,"mode":"ojs"},{"id":381,"value":"R0ByInfectorAge = {\n let R0ByAge = { } \n infectors.forEach((d,i) => \n R0ByAge[d.age] = (R0ByAge[d.age] || 0) + infectionStats[i].cumulativeInfections\n )\n Object.keys(R0ByAge).forEach(age => \n R0ByAge[age] = R0ByAge[age] / infectors.filter(d => d.age == age).length\n )\n return R0ByAge\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":380,"value":"md `R0 by infector symptomaticity is cumulative infections for each infector whether they are asymptomatic / total number of infectors that have the same symptomaticity. Note that this does not adress infections from presymptomatic individuals.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":379,"value":"R0ByInfectorSymptomaticity = {\n let R0BySymptomaticity = { } \n infectors.forEach((d,i) => \n R0BySymptomaticity[d.asymptomatic ? 'Asympt':'Symptomatic'] = (R0BySymptomaticity[d.asymptomatic ? 'Asympt':'Symptomatic'] || 0) + infectionStats[i].cumulativeInfections\n )\n Object.keys(R0BySymptomaticity).forEach(symptomaticity => \n R0BySymptomaticity[symptomaticity] = R0BySymptomaticity[symptomaticity] / infectors.filter(d => symptomaticity == (d.asymptomatic ? 'Asympt':'Symptomatic') ).length\n )\n \n return R0BySymptomaticity\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":67,"value":"md`## Caveats \n- I am not an epidemiologist or infectious disease expert. While I believe these calculations correctly translate our game simulation parameters to an average R0, I may well have errors. \n- R0 is an important parameter, but how quickly R changes is also very important, and not covered here.\n- This notebook was put together under significant time-constraints. There's plenty more cleaning up that should happen, but figured it was more useful to put out even with additional enhancements listed below`","pinned":false,"hidden":false,"mode":"ojs"},{"id":343,"value":"md`## Potential Enhancements \n- Add breakdowns to show number of infections by person per age groups, during presymptomatic period, √ for people who were asymptomatic vs symptomatic people, √ to housemates vs out in public, \n- Add plots to show distribution of infectious caused instead of just mean\n- Allow infectiousness to be adjustable over the over disease course\n- For help with performance don't run simulation until parameters are all changed and recalc of R0 is triggered with a button or something. \n- Make parameters more visually compact and easier to adjust to see most things on one screen\n- Adjust number of destinations per day for whole population (which both increases the number of times infectious people go out and the number of other people they find at those places)`","pinned":false,"hidden":false,"mode":"ojs"},{"id":328,"value":"md`## Appendix`","pinned":false,"hidden":false,"mode":"ojs"},{"id":248,"value":"import {slider,select} from \"@jashkenas/inputs\"","pinned":true,"hidden":false,"mode":"ojs"},{"id":104,"value":"d3 = require(\"d3@5\")","pinned":true,"hidden":false,"mode":"ojs"},{"id":376,"value":"_ = require(\"lodash\")","pinned":true,"hidden":false,"mode":"ojs"},{"id":375,"value":"roundProb = function(x) {return Math.floor(x) + (Math.random() < (x - Math.floor(x))) } ","pinned":true,"hidden":false,"mode":"ojs"},{"id":374,"value":"consoleLogging = 0 // if this changes to 1 there's a whole lot of console-logging in the infectionsByInfector cell -- which can get slow when 1K infectors. When I wanted to inspect the console log I'd set numInfectors to 7-10 first. ","pinned":true,"hidden":false,"mode":"ojs"}],"files":"$W34","stdlib":"1","titling":"auto"}