R version 2.8.0 (2008-10-20)
Copyright (C) 2008 The R Foundation for Statistical Computing
ISBN 3-900051-07-0

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

  Natural language support but running in an English locale

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> par4 = '0.6'
> par3 = '15'
> par2 = '45'
> par1 = '3650'
> par1 <- as.numeric(par1)
> par2 <- as.numeric(par2)
> par3 <- as.numeric(par3)
> par4 <- as.numeric(par4)
> numsuccessbig <- 0
> numsuccesssmall <- 0
> bighospital <- array(NA,dim=c(par1,par2))
> smallhospital <- array(NA,dim=c(par1,par3))
> bigprob <- array(NA,dim=par1)
> smallprob <- array(NA,dim=par1)
> for (i in 1:par1) {
+ bighospital[i,] <- sample(c('F','M'),par2,replace=TRUE)
+ if (as.matrix(table(bighospital[i,]))[2] > par4*par2) numsuccessbig = numsuccessbig + 1
+ bigprob[i] <- numsuccessbig/i
+ smallhospital[i,] <- sample(c('F','M'),par3,replace=TRUE)
+ if (as.matrix(table(smallhospital[i,]))[2] > par4*par3) numsuccesssmall = numsuccesssmall + 1
+ smallprob[i] <- numsuccesssmall/i
+ }
> tbig <- as.matrix(table(bighospital))
> tsmall <- as.matrix(table(smallhospital))
> tbig
   [,1]
F 82300
M 81950
> tsmall
   [,1]
F 27262
M 27488
> numsuccessbig/par1
[1] 0.06027397
> bigprob[par1]
[1] 0.06027397
> numsuccesssmall/par1
[1] 0.1476712
> smallprob[par1]
[1] 0.1476712
> numsuccessbig/par1*365
[1] 22
> bigprob[par1]*365
[1] 22
> numsuccesssmall/par1*365
[1] 53.9
> smallprob[par1]*365
[1] 53.9
> postscript(file="/var/www/html/freestat/rcomp/tmp/1r8xq1286814547.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) 
> plot(bigprob,col=2,main='Probability in Large Hospital',xlab='#simulated days',ylab='probability')
> dev.off()
null device 
          1 
> postscript(file="/var/www/html/freestat/rcomp/tmp/2khxs1286814547.ps",horizontal=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) 
> plot(smallprob,col=2,main='Probability in Small Hospital',xlab='#simulated days',ylab='probability')
> dev.off()
null device 
          1 
> 
> #Note: the /var/www/html/freestat/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab
> load(file="/var/www/html/freestat/rcomp/createtable")
> 
> a<-table.start()
> a<-table.row.start(a)
> a<-table.element(a,'Exercise 1.13 p. 14 (Introduction to Probability, 2nd ed.)',2,TRUE)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'Number of simulated days',header=TRUE)
> a<-table.element(a,par1)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'Expected number of births in Large Hospital',header=TRUE)
> a<-table.element(a,par2)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'Expected number of births in Small Hospital',header=TRUE)
> a<-table.element(a,par3)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'Percentage of Male births per day<br />(for which the probability is computed)',header=TRUE)
> a<-table.element(a,par4)
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'#Females births in Large Hospital',header=TRUE)
> a<-table.element(a,tbig[1])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'#Males births in Large Hospital',header=TRUE)
> a<-table.element(a,tbig[2])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'#Female births in Small Hospital',header=TRUE)
> a<-table.element(a,tsmall[1])
> a<-table.row.end(a)
> a<-table.row.start(a)
> a<-table.element(a,'#Male births in Small Hospital',header=TRUE)
> a<-table.element(a,tsmall[2])
> a<-table.row.end(a)
> a<-table.row.start(a)
> dum1 <- paste('Probability of more than', par4*100, sep=' ')
> dum <- paste(dum1, '% of male births in Large Hospital', sep=' ')
> a<-table.element(a, dum, header=TRUE)
> a<-table.element(a, bigprob[par1])
> a<-table.row.end(a)
> dum <- paste(dum1, '% of male births in Small Hospital', sep=' ')
> a<-table.element(a, dum, header=TRUE)
> a<-table.element(a, smallprob[par1])
> a<-table.row.end(a)
> a<-table.row.start(a)
> dum1 <- paste('#Days per Year when more than', par4*100, sep=' ')
> dum <- paste(dum1, '% of male births occur in Large Hospital', sep=' ')
> a<-table.element(a, dum, header=TRUE)
> a<-table.element(a, bigprob[par1]*365)
> a<-table.row.end(a)
> dum <- paste(dum1, '% of male births occur in Small Hospital', sep=' ')
> a<-table.element(a, dum, header=TRUE)
> a<-table.element(a, smallprob[par1]*365)
> a<-table.row.end(a)
> a<-table.end(a)
> table.save(a,file="/var/www/html/freestat/rcomp/tmp/3jrtp1286814547.tab") 
> try(system("convert tmp/1r8xq1286814547.ps tmp/1r8xq1286814547.png",intern=TRUE))
character(0)
> try(system("convert tmp/2khxs1286814547.ps tmp/2khxs1286814547.png",intern=TRUE))
character(0)
> 
> 
> proc.time()
   user  system elapsed 
  6.067   1.274   6.168