I started working with R recently and need to do a simulation for 100 000 people that divides them in exposition factor for a disease, gender and presence or absence of disease. My background info is that the probability of being exposed to this risk factor is 0.15 in this population. The probability of being a woman, among exposed people is 0.2 and among non exposed people is 0.5. I also know that there are 4 gamma distributions - divided by gender and exposition factor - that characterize the disease. If X> 35, from each gamma, then you have the disease. I used the following code to answer this question:

#vector that represents non exposed ("non exposed" = 1, "exposed" = 0)
gender<-rep(0,n) #vector that reprsents gender ("women" = 1, "men" = 0)
#vector that represents disease ("with disease" = 0, "without disease" = 1)

#in 100 000 cases atributes non expostion
for(i in 1:n){if(u[i]>0.15)exposicao[i]<-1}

#cycle that determines the gender of the subjects
for(j in 1:n){{if(exposicao[j]==1 & b[j]>0.5)genero[j]<-1} 
{if(exposicao[j]==0 & b[j]<0.2)genero[j]<-1}} 

#Cycle that determines the state of disease from the subjects. I define people with
disease as people with c[x] inferior to the probability of having the disease
(P(X>35)), according to their gamma distribution. I used the results from that

     for(x in 1:n){{if(expositon[x]==1 & gender[x]==1 & c[x]<0.06840084)doentes[x]<-1} 
 #if the person is non exposed and woman X2~G(12,2)
 {if(exposition[x]==0 & gender[x]==1 & c[x]<0.1698673)doentes[x]<-1} 
 #if it's exposed and woman X3~G(14,2)
 {if(exposicao[x]==1 & genero[x]==0 & c[x]<0.02010428)doentes[x]<-1} 
 #if it's non exposed and man X1~G(10,2)
 {if(exposicao[x]==0 & genero[x]==0 & c[x]<0.3275424)doentes[x]<-1}}
 #if it's exposed and man X4~G(16,2)

1 Answer 1


There are a few improvements to your (I'm assuming that genero and gender are the same variable).

  1. Very important. Use a space around the assignment operator! x<-5 Does this mean assign 5 to x to is x less than -5
  2. Rather then rep(0,n) just use numeric(n)
  3. Your loops should be rewritten. For example:

    ##I presume that exposicao is exposition.
    for(i in 1:n){if(u[i]>0.15)exposicao[i]<-1}

    can be written as:

    exposition[runif(n) > 0.15] = 1


    for(j in 1:n){{if(exposicao[j]==1 & b[j]>0.5)genero[j]<-1} 
    {if(exposicao[j]==0 & b[j]<0.2)genero[j]<-1}} 

    can be written as:

    genero[exposicao==1 & b>0.5] = 1
    genero[exposicao==0 & b>0.2] = 1

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.