Simulation — Rodent Plague in Manhattan
A simulation of rodent plague influence and the emergency response in Manhattan
COURSE: Datamining the City — Session A| GSAPP Fall 2018 | Columbia University
TEAM: Yingjun Mou, Shuyang Huang, Yue Dong
We aim to do a project simulating the spread of potential rodent plague in Manhattan. With 311 complaint data, we are able to recognize the locations where the sources of infection might be. The population dynamics simulation in Processing will be the tool for us to illustrate the spreading trend of the epidemics and help us to know how to apply direct and immediate response toward the emergent accident. This report aims to address the following research questions:
- Who will ultimately get infected according to the existing condition and the assumption?
- What is the minimum number of infected people that will result in the whole neighborhood getting infected?
- How will the cure rate and human fluidity influence the spread of the plague?
Regarding the rodent plague as a case study, this prototype is to serve as a template for other epidemics.
Background:
A research in 2014 shows that there are around 2 million rats in total in New York City, which is approximately 24% of the number of humans.
Furthermore, rodent problem is also one of the top concerns for New Yorkers. From 2017 to current, around 876,000 complaints are recorded through the 311 system in Manhattan. Rodent complaints, including rat or mice sightings and signs of rodents, is the Top 17 concern and comprise around 1.6% of total records.
So what if a rodent plague outbreak?
Data Input:
A. Demographic Data of Manhattan
B. The modeling of the city
C. History record of rodent presence
Data Output:
A. Calculate the “Critical area”. i.e. if there is limited medical resources (imagine a new type of panacea which can cure infected people instantly), in which neighborhood should government use the panacea, in order to stop the spreading of disease, and minimize the number of infected people?
B. Visualize the general distribution of infected people throughout Manhattan.
C. Show whether the given population fluidity and cure rate can curb the spreading of disease or aggravate it
Data Structure:
We will have two interfaces: City and People. And while the City has three sub-classes, the class People will have three different states.
A. City:
1. Neighborhood — — medium population number, low dynamic level
2. Street — — high population number, high dynamic level
3. Park — — low population number, low dynamic level
B. People:
There are 4 different states of a person (represented by the color of pixel)
1. Healthy (Green pixel)
2. Being infected (Yellow pixel, and will become red the next day)
3. Infected (Red pixel)
4. Dead (Pixel removed)
Rule of Infection:
(1) If there are >= 2 out of 4 nearby people infected, then the persons will be infected the next day (each people is represented as a square pixel)
(2) After 50 days (one day each frame), an infected person will die
(3) There is a preset probability of being cured, and a preset probability of moving around(and to spread the disease).
Methodology:
(1) Data preparation:
We use R studio and GIS to clean and aggregate the data. 311 rodent complaint dataset on the New York open data is made by Department of Health and Mental Hygiene (DOHMH) and is updated daily since 12/16/2015. To simplify the analyzing process, we extract the complaint data from 2017 to 2018 and spatial join them into census blocks. And several areas are recognized to have more rodent inspection, such as Southern Bronx, Upper West Side, East Village, and Lower East Side.
(2) Simulation:
Manhattan has been divided into 27x72 2D grids. CSV made in previous step is imported in processing and read line-by-line, providing each grid with attribute of its own population density and plague outbreak probability. In each grid, 25 pixels are used to represent the population in each block referring to the ACS data. And locations with higher population density and rodent records have been imported in the processing canvas and acted as the source of infection. Disease cure rate and population fluidity are regarded as the key factors of terminating the spread of plague.
Disease cure rate in this report is presented as percentage format. The 2.5% disease cure rate indicates that in 100 people randomly 2.5 people can be cured in an iteration.
Population fluidity is also presented as percentage format. The 50% population fluidity demonstrates that in 100 people, 50 people would randomly move in an iteration.
Different disease cure rates and population fluidity are set in the research process to test the best solution toward quenching the spread of plague. We define 2.5% cure rate as a low level cure rate, 25% cure rate as a high level cure rate; and 10% population fluidity as a low population fluidity, 50% population fluidity as a high population fluidity.
Results:
Below are 5 different scenarios that were simulated based on different combinations of given disease cure rates and population fluidity.
Conclusion:
(1) Based on the history record of rodent, and the population distribution in Manhattan, there are several areas which has a higher probability of suffering from rodent plague, such as Southern Bronx, Upper West Side, Midtown West, East Village, and Lower East Side.
(2) In terms of quenching the spread of plague, the Higher Cure Rateworks better than the Lower Population Fluidity. Population Fluidity is more about how fast the plague will spread.
(3) Assume that the Population Fluidity is 10%(means there are 10% of the total population that constantly moving around), the city need to achieve a cure rate of at least 10% to stop the plague from spreading. Assume that the Population Fluidity is 50%, then a cure rate of 25% is enough for quenching the plague in the city.
Data Processing steps:
(1) Data cleaning code
library(readr)
library(lubridate)
Rodent <- read_csv("data/Rodent_Inspection.csv")
summary(Rodent)
library(dplyr)
Rodent %>% select(`ZIP_CODE`, `X_COORD`, `Y_COORD`, `LATITUDE`,`LONGITUDE`, `BOROUGH`,
`INSPECTION_DATE`, `LOCATION`) %>%
filter(complete.cases(.)) %>% filter(`ZIP_CODE` != 0 & `X_COORD`!= 0 &
`Y_COORD`!= 0 & `LATITUDE`!= 0 & `LONGITUDE` != 0)-> ComplaintRodent
names(ComplaintRodent) <- c("Zip", "XCoord", "YCoord", "Lat", "Lon", "Borough", "Date", "Location")
View(ComplaintRodent)library(stringr)
ComplaintRodent %>% mutate (Date2 = str_sub(ComplaintRodent$Date, 1, 19)) %>%
mutate(Date3 = as.POSIXct(ComplaintRodent$Date, format= "%m/%d/%Y %H:%M:%S"))%>%
mutate(ymd = paste0(str_pad(year(Date3), 4, "left", "0"),
str_pad(month(Date3), 2, "left", "0"),
str_pad(day(Date3), 2, "left", "0"))) %>%
select(ymd, Location, Zip, Borough, XCoord, YCoord, Lat, Lon, Date3)-> ComplaintRodent2ComplaintRodent3 <- subset(ComplaintRodent2, Date3 >= "2017-01-01")
write.csv(ComplaintRodent3, file = "Rodent2017.csv",row.names=FALSE)
(2) Data Mining and Simulation:
from random import randint
import csv
import operator
import collections
class Person(object):
#=========CONSTRUCTORS==========
#x if the row index, y is the column index, state is an int indicate the health condition
#state=0 healthy(green), state=1 being infected(yellow), state=2 infected state=2+50=52 dead
#Currently the time between getting infected to dead is set to 50===========================
def __init__(self, state):
self.state = state
#=========OBSERVERS==========
#return the health state of a person
def checkState(self):
return self.state
#=========MODIFIERS==========
#return the health state of a person
def changeState(self, my_state):
self.state = my_state
class Population(object):
#=========CONSTRUCTORS==========
def __init__(self, maxX, maxY):
self.maxX = maxX
self.maxY = maxY
row = [None]*maxX
self.data = []
for i in range(maxX):
row = []
for j in range(maxY):
row.append([None])
self.data.append(row)
#=========OBSERVERS==========
#count the number of adjacent infected people
def countAdj(self, x, y):
count =0
#case 0: no person at x,y
if (self.data[x][y]==None):
return -1
#case 1: x==0, y==0
elif (x==0 and y==0):
if (self.data[x+1][y] != [None] and self.data[x+1][y].checkState() >= 2):
count += 1
if (self.data[x][y+1] != [None] and self.data[x][y+1].checkState() >= 2):
count += 1
return count
#case 2: x==maxX, y==maxY
elif (x==self.maxX-1 and y==self.maxY-1):
if (self.data[x-1][y] != [None] and self.data[x-1][y].checkState() >= 2):
count += 1
if (self.data[x][y-1] != [None] and self.data[x][y-1].checkState() >= 2):
count += 1
return count
#case 3: x==0, y==maxY
elif (x==0 and y==self.maxY-1):
if (self.data[x+1][y] != [None] and self.data[x+1][y].checkState() >= 2):
count += 1
if (self.data[x][y-1] != [None] and self.data[x][y-1].checkState() >= 2):
count += 1
return count
#case 4: x==maxX, y==0
elif (x==self.maxX-1 and y==0):
if (self.data[x-1][y] != [None] and self.data[x-1][y].checkState() >= 2):
count += 1
if (self.data[x][y+1] != [None] and self.data[x][y+1].checkState() >= 2):
count += 1
return count
#case 5: x==0
elif (x==0):
if (self.data[x+1][y] != [None] and self.data[x+1][y].checkState() >= 2):
count += 1
if (self.data[x][y-1] != [None] and self.data[x][y-1].checkState() >= 2):
count += 1
if (self.data[x][y+1] != [None] and self.data[x][y+1].checkState() >= 2):
count += 1
return count
#case 6: x==maxX
elif (x==self.maxX-1):
if (self.data[x-1][y] != [None] and self.data[x-1][y].checkState() >= 2):
count += 1
if (self.data[x][y-1] != [None] and self.data[x][y-1].checkState() >= 2):
count += 1
if (self.data[x][y+1] != [None] and self.data[x][y+1].checkState() >= 2):
count += 1
return count
#case 7: y==maxY
elif (y==self.maxY-1):
if (self.data[x-1][y] != [None] and self.data[x-1][y].checkState() >= 2):
count += 1
if (self.data[x][y-1] != [None] and self.data[x][y-1].checkState() >= 2):
count += 1
if (self.data[x+1][y] != [None] and self.data[x+1][y].checkState() >= 2):
count += 1
return count
#case 8: y==0
elif (y==0):
if (self.data[x-1][y] != [None] and self.data[x-1][y].checkState() >= 2):
count += 1
if (self.data[x][y+1] != [None] and self.data[x][y+1].checkState() >= 2):
count += 1
if (self.data[x+1][y] != [None] and self.data[x+1][y].checkState() >= 2):
count += 1
return count
#case 9: others
else:
if (self.data[x-1][y] != [None] and self.data[x-1][y].checkState() >= 2):
count += 1
if (self.data[x+1][y] != [None] and self.data[x+1][y].checkState() >= 2):
count += 1
if (self.data[x][y-1] != [None] and self.data[x][y-1].checkState() >= 2):
count += 1
if (self.data[x][y+1] != [None] and self.data[x][y+1].checkState() >= 2):
count += 1
return count
#=========MODIFIERS==========
def infect(self,x,y):
#Pre-condition: data[x][y] != None
if(self.data[x][y] != [None]):
#case 1: from healthy to being infected
if (self.data[x][y].checkState()==0 and self.countAdj(x,y)>=2):
self.data[x][y].changeState(1)
#case 2: from being infected to infected
elif (self.data[x][y].checkState()==1):
self.data[x][y].changeState(2)
#case 3: after getting infected, disease becoming aggrevated if it's not cured
elif (self.data[x][y].checkState()>=2 and self.data[x][y].checkState()<52):
self.data[x][y].changeState(self.data[x][y].checkState()+1)
else:
return
def cured(self,x,y):
if(self.data[x][y] != [None]):
prob = randint(0,4)
if(prob==0):
if (self.data[x][y].checkState()>0): # and self.data[x][y].checkState()<52
self.data[x][y].changeState(0)def move(self, x, y):
#randomly swap population pixels,
if(density_map[x/5][y/5] != 0.0 and self.data[x][y] == [None]):
direction = [0,1,2,3]
if(y-1>=0):
if(self.data[x][y-1] == [None] or density_map[x/5][(y-1)/5] <0.3):
direction.remove(0)
if(y+1<=self.maxY):
if(self.data[x][y+1] == [None] or density_map[x/5][(y+1)/5] <0.3):
direction.remove(1)
if(x-1>=0):
if(self.data[x-1][y] == [None] or density_map[(x-1)/5][y/5] <0.3):
direction.remove(2)
if(x+1<self.maxX):
if(self.data[x+1][y] == [None] or density_map[(x+1)/5][y/5] <0.3):
direction.remove(3)
if(x+1>=self.maxX):
direction.remove(3)prob = randint(0,2)
if(len(direction) >0 and prob<1):
l = len(direction)
random_index = randint(0,l-1)
random_direction = direction[random_index]
if(random_direction==0):
self.data[x][y] = self.data[x][y-1]
self.data[x][y-1] = [None]
elif(random_direction==1):
self.data[x][y] = self.data[x][y+1]
self.data[x][y+1] = [None]
elif(random_direction==2):
self.data[x][y] = self.data[x-1][y]
self.data[x-1][y] = [None]
else:
if(x+1<self.maxX):
self.data[x][y] = self.data[x+1][y]
self.data[x+1][y] = [None]
#INITIALIZATION OF THE POPULATION SYSTEM
global cell_width, cell_gap, my_population, denstiy_map
canvas_width = 1600
canvas_height = 600
cell_width = 4
cell_gap = 1
density_map = []#Read CSV file and load the density data
with open('density and outbreak.csv') as csv_file:
csv_reader = csv.reader(csv_file, delimiter=',')
line_count = 0
current_i = -1
for row in csv_reader:
if line_count != 0:
if int(row[0]) != current_i:
density_map.append([])
current_i = int(row[0])
density_map[current_i].append(float(row[2]))
line_count += 1
unit_size = cell_width + cell_gap
num1 = canvas_width / unit_size
num2 = canvas_height / unit_size
my_population = Population(num1, num2)
for i in range(num1):
for j in range(num2):
#Probability of the outbreak of plauge
if(density_map[i/5][j/5] != 0):
threshold = int(float(density_map[i/5][j/5])*10)
prob = randint(0, 10)
if(prob <= threshold):
outbreak = randint(0,10)
if(outbreak==0):
my_population.data[i][j] = Person(2)
else:
my_population.data[i][j] = Person(0)def setup():
frameRate(96)
size (canvas_width,canvas_height)
smooth()
background(0)def draw():
global cell_width, cell_gap, my_population,img
stroke (127)
strokeWeight(1)
strokeCap(SQUARE)
background(0)
#display of the population
for i in range(num1):
for j in range(num2):
#CHECK NULL
if (my_population.data[i][j] != None and my_population.data[i][j] != [None]):
if (my_population.data[i][j].checkState() == 0):
fill(0,255,0)
rect(i*unit_size, j*unit_size, cell_width, cell_width)
elif(my_population.data[i][j].checkState() == 1):
fill(255,255,0)
rect(i*unit_size, j*unit_size, cell_width, cell_width)
#ADD DIFFERENT GRADIENT OF RED
elif(my_population.data[i][j].checkState() >= 2 and my_population.data[i][j].checkState() < 52):
delta = my_population.data[i][j].checkState()-2
fill(255-delta*5,0,0)
rect(i*unit_size, j*unit_size, cell_width, cell_width)
elif(my_population.data[i][j].checkState() ==52):
fill(127)
rect(i*unit_size, j*unit_size, cell_width, cell_width)
else:
print("Exception")
for i in range(num1):
for j in range(num2):
my_population.infect(i,j)
my_population.cured(i,j)
for i in range(num1):
for j in range(num2):
my_population.move(i,j)
title()
saveFrame("outputE/plagueE_####.png")
def title():
fill(255);
textSize(16);
text("Cure Rate=25% \nPopulation Fluidity=50%", 30, 30);
