Skip to content

Instantly share code, notes, and snippets.

@rcsmit
Last active August 15, 2022 22:33
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save rcsmit/4c4e32f738855c159a0e61e17a4e449b to your computer and use it in GitHub Desktop.
Save rcsmit/4c4e32f738855c159a0e61e17a4e449b to your computer and use it in GitHub Desktop.
disabled_with_various_infected
for (m in 1:10){
x = m/10 #new in R so I dont know how to make steps of .1 in loops ;)
monsterlist <- c() #This will hold the number of disabled in each sim run
for (h in 1:1){ #Repeat simulation 1 times
disabled <- c() #This will hold the disabled within each sim run
population <- as.vector(c(rep(0,10000))) #The population of 10,000
for(j in 1:30) { # Number of Years
infected <- sample(1:length(population),floor(x*length(population))) # infect x% of the population every year
population[infected] <- population[infected]-1 #infected get an increasingly-negative counter
# The next 2 paragraphs control who gets LC. One should be commented out.
# The first puts everyone in the same risk category.
# The second creates low, medium, and high risk categories.
for (k in 1:length(population)){ #Use this paragraph if risk to everyone of LC is equal.
#Use next paragraph if LC is based on low risk, medium risk, and high risk groups.
# #if (runif(1) <= pweibull(abs(population[k]),scale=6,shape=2)){ #Use this one if successive Covid infections increase LC risk
if (runif(1) <= 0.1) { #Use this one if risk of LC is constant. Comment out one of these.
population[k] <- population[k]+ 100 #positive number means Long Covid
}
}
# for (p in 1:length(population)){
# if (p <= 3000) { #The first 3000 in the population are "low risk" for LC.
# if(runif(1) <= 0.05) { #Only a 5% chance of getting LC with each infection.
# population[p] <- population[p] + 100
# }
# } else if (p >= 7000) { #High risk group for LC.
# if(runif(1) <= 0.25) { #This group has a 25% chance of getting LC with each infection.
# population[p] <- population[p] + 100
# }
# } else { #This is the "medium" risk group of getting LC with each infection.
# if(runif(1) <= 0.1) { #This group has a 10% chance of getting LC with each infection.
# population[p] <- population[p] + 100
# }
# }
# }
for (n in 1:length(population)){ # give them a 50% chance to recover from LC this year
if (population[n] > 0){
if(runif(1)>=0.5){
population[n] <- population[n] - 100
}
}
}
del_list <- c()
for (m in 1:length(population)){
if (population[m] > 0){
if(runif(1)>=0.9){
disabled <- c(disabled,population[m])
del_list <- c(del_list,m)
}
}
}
if (length(del_list)>0){
population <- population[-del_list]
}
}
monsterlist <- c(monsterlist,length(disabled))
print (paste("infected:", r, " result" , sum(monsterlist)/1000, "%")) #gives the % (not the decimal) average over 100 simulations
}
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment