worldCurrencies
worldCurrencies

Reputation: 467

How to get this probability outcome in R?

New to R. Not sure how to go about this problem. Essentially I want to create a Monte Carlo simulation. However this is how it's supposed to go: There are only 3 people (A,B,C) And only come into contact once in this order: (A-->B) (B-->C). Starting out, A is 100% sick, while B and C are 0% sick. Each time a person comes into contact with another, they have a 50% chance of being sick How would i go and replicate that in R? I understand the mathematics behind this but unsure how to code it. This is what i have so far:

a='positive'
corona = c('positive','negative')
sample(x = corona, size = 1, replace = TRUE)

I know the output will only give me the results of whether B is sick or not. How would i continue to see if C is sick?

Upvotes: 1

Views: 193

Answers (1)

MrFlick
MrFlick

Reputation: 206197

You could write a helper function to do the "infecting". Let's assume that 1 is "sick" and 0 is "not sick"

infect <- function(x, rate=.5) {
   current_sick <- x==1
   new_sick <- runif(sum(current_sick)) < rate
   x[current_sick] <- as.numeric(new_sick) # turn TRUE/FALSE to 1/0
   x
}

This function looks where all the 1's are, then for each of those individuals, it tosses a coin to see of the next person will be sick (not changing any of the non-sick values)

Then to test with 1000 people, you can do

A <- rep(1, 1000) # all sick
B <- infect(A)
C <- infect(B)

This draws 1000 A's. And we go though two transmission steps, doing all 1000 samples at once. If you run mean(C) you should get something close to 0.25.

Upvotes: 3

Related Questions