set.seed(984293781) media=12 dp=1 N=5000 hemo=rnorm(N,media,dp) #nivel de hemoglobina em mulheres jovens e saud?veis #jpeg(file='taxashemog.jpg') hist(hemo,main='',xlab='Taxas de hemoglobina',freq=FALSE) #dev.off() N=length(hemo) media.pop=round(mean(hemo)) dp.pop=round(sd(hemo)) #Sobrepor a distribui??o de probabilidade normal d=dnorm(seq(min(hemo),max(hemo),l=200),media.pop,dp.pop) lines(seq(min(hemo),max(hemo),l=200),d,col=2) # Na pr?tica a m?dia e o desvio-padr?o s?o desconhecidos!!! # # Ent?o uma amostra aleat?ria ? selecionada e com base nesta amostra # temos uma estimativa da m?dia e o desvio-padr?o. # #Selecionando aleatoriamente amostra de tamanho n n=10 grupo=sample(1:N,size=n,replace=TRUE) amostra=hemo[grupo] round(amostra,2) round(mean(amostra),2) round(sd(amostra),2) #Selecionando outras n mulheres...temos um resultado diferente... grupo=sample(1:N,size=n,replace=TRUE) amostra=hemo[grupo] round(amostra,2) round(mean(amostra),2) round(sd(amostra),2) # PERGUNTAS: # COMO ESTIMAR A M?DIA E CALCULAR A PRECIS?O DA ESTIMATIVA? # SER? QUE EXISTE UM PADR?O ESPERADO NO COMPORTAMENTO DAS M?DIAS? # Vamos tentar responder estas perguntas com um exerc?cio de simula??o. # Selecionando repetidamente grupos de n mulheres # e calculando as m?dias e desvios-padr?o n=25 nsim=1000 medias=NULL dps=NULL for (i in 1:nsim){ grupo=sample(1:N,size=n,replace=TRUE) amostra=hemo[grupo] medias=c(medias,mean(amostra)) dps=c(dps,sd(amostra)) } par(mfrow=c(1,2)) hist(hemo,main='Histograma do n?vel de hemoglobina',xlim=c(min(hemo),max(hemo)),freq=FALSE) hist(medias,main='Histograma das m?dias (n=10)',xlim=c(min(hemo),max(hemo)),freq=FALSE) mean(medias) #m?dia das m?dias media.pop #m?dia dos n?veis de hemoglobina sd(medias) #desvio-padr?o das m?dias dp.pop #desvio-padr?o dos n?veis de hemoglobina dp.pop/sqrt(n) # Note que o resultado desta conta se aproxima do sd(medias) # ou seja, o desvio-padr?o das m?dias ? dado pelo # desvio-padr?o da popula??o dividido por sqrt(n). # Resultado te?rico: media.amostral~N(media.pop,dp.pop/sqrt(n)) par(mfrow=c(1,1)) dd=dnorm(seq(min(medias),max(medias),l=200),media.pop,dp.pop/sqrt(n)) hist(medias,freq=FALSE,xlab='M?dias amostrais') lines(seq(min(medias),max(medias),l=200),dd,col=2) # Com este resultado ? poss?vel obter uma estimativa intervalar # para a m?dia populacional. # demonstrar no quadro #V?rios intervalos de confian?a num gr?fico nsim=100 li=medias-1.96*dps/sqrt(n) ls=medias+1.96*dps/sqrt(n) plot(1:nsim,medias[1:nsim],ylim=c(min(li[1:nsim]),max(ls[1:nsim])),xlab='Amostra',ylab='Intervalo de Confian?a') abline(h=media.pop,lty=3) segments(1:nsim,li[1:nsim],1:nsim,ls[1:nsim]) ########################## #TEOREMA CENTRAL DO LIMITE ########################## dat=c(9,1,0,7,5,6,9,5,8,8,1,0,5,7,6,5,0,2,1,2,1,8,8,8,5,2,4,8,3,1,6,5,5,7,4,1,7,3,3,3,2,8,1,8,5,8,4,0,1,9,2,1,6,9,4,4,7,6,1,7,1,9,7,9,7,2,7,7,0,8,1,6,3,8,0,5,7,4,8,6,7,0,2,8,8,7,2,5,4,1,8,6,8,3,5,8,2,7,2,4) hist(dat,breaks=seq(0,10,by=1),include.lowest = TRUE, right = FALSE,main='',xlab='x') medias=NULL for (i in 1:200){ x=dat[sample(1:length(dat),size=40,replace=TRUE)] medias=c(medias,mean(x)) } medias hist(medias,freq=FALSE,main='',xlab='m?dias') #Selecionando aleatoriamente amostra de tamanho 6 n=6 grupo=sample(1:N,size=n,replace=TRUE) amostra=hemo[grupo] round(amostra,2) round(mean(amostra),2) round(sd(amostra),2) #Selecionando outras 6 mulheres...temos um resultado diferente... grupo=sample(1:N,size=n,replace=TRUE) amostra=hemo[grupo] round(amostra,2) round(mean(amostra),2) round(sd(amostra),2) # PERGUNTAS: # COMO ESTIMAR A M?DIA E CALCULAR A PRECIS?O DA ESTIMATIVA? # SER? QUE EXISTE UM PADR?O ESPERADO NO COMPORTAMENTO DAS M?DIAS? # Vamos tentar responder estas perguntas com um exerc?cio de simula??o. #Selecionando repetidamente grupos de 6 mulheres #e calculando as m?dias e desvios-padr?o n=5 nsim=1000 medias=NULL dps=NULL for (i in 1:nsim){ grupo=sample(1:N,size=n,replace=TRUE) amostra=hemo[grupo] if(i<=10) print(round(amostra,2)) medias=c(medias,mean(amostra)) dps=c(dps,sd(amostra)) } par(mfrow=c(1,2)) #jpeg(file='histhemo.jpg') hist(hemo,main='Histograma do n?vel de hemoglobina',xlim=c(min(hemo),max(hemo)),freq=FALSE) #dev.off() #jpeg(file='histmedias6.jpg') hist(medias,main='Histograma das m?dias (n=6)',xlim=c(min(hemo),max(hemo)),freq=FALSE) #dev.off() mean(medias) #m?dia das m?dias media.pop #m?dia dos n?veis de hemoglobina sd(medias) #desvio-padr?o das m?dias dp.pop #desvio-padr?o dos n?veis de hemoglobina dp.pop/sqrt(n) # Note que o resultado desta conta se aproxima do sd(medias) # ou seja, o desvio-padr?o das m?dias ? dado pelo # desvio-padr?o da popula??o dividido por sqrt(n). # Resultado te?rico: media.amostral~N(media.pop,dp.pop/sqrt(n)) # Z=(media.amostral-media.pop)/(dp.pop/sqrt(n))~N(0,1) par(mfrow=c(1,1)) z=(medias-media.pop)/(dp.pop/sqrt(n)) dd=dnorm(seq(-5,5,l=200),0,1) hist(z,freq=FALSE,xlim=c(-5,5)) lines(seq(-5,5,l=200),dd,col=2) # Com este resultado ? poss?vel obter uma estimativa intervalar # para a m?dia populacional. # Na pr?tica dp.pop tamb?m n?o ? conhecido... # ...ent?o se estimarmos dp.pop usando o desvio-padr?o da amostra # t=(media.amostral-media.pop)/(dp.amostral/sqrt(n)) ~ t-Student t=(medias-media.pop)/(dps/sqrt(n)) tt=dt(seq(-5,5,l=200),n-1) hist(t,freq=FALSE,xlim=c(-5,5)) lines(seq(-5,5,l=200),dd,col=2) lines(seq(-5,5,l=200),tt,col=4) ########################################################################## # Voltando ao exemplo do n?vel de hemoglobina em mulheres jovens saud?veis #Selecionando aleatoriamente amostra de tamanho 10 n=10 grupo=sample(1:N,size=n,replace=TRUE) amostra=hemo[grupo] # Com base nesta amostra, o intervalo de confian?a de 95% # para o n?vel m?dio de hemoglobina ser?: li=mean(amostra)-qt(0.975,n-1)*sd(amostra)/sqrt(n) ls=mean(amostra)+qt(0.975,n-1)*sd(amostra)/sqrt(n) cbind(li,ls) amplitude=ls-li amplitude # Se usarmos a distribui??o normal: # o intervalo fica mais estreito mas n?o ? o mais correto li=mean(amostra)-qnorm(0.975,0,1)*sd(amostra)/sqrt(n) ls=mean(amostra)+qnorm(0.975,0,1)*sd(amostra)/sqrt(n) cbind(li,ls) amplitude=ls-li amplitude #V?rios intervalos de confian?a num gr?fico nsim=100 medias=NULL dps=NULL for (i in 1:nsim){ grupo=sample(1:N,size=n,replace=TRUE) amostra=hemo[grupo] if(i<=10) print(round(amostra,2)) medias=c(medias,mean(amostra)) dps=c(dps,sd(amostra)) } li=medias-qt(0.975,n-1)*dps/sqrt(n) ls=medias+qt(0.975,n-1)*dps/sqrt(n) plot(1:nsim,medias,ylim=c(min(li),max(ls)),xlab='Amostra',ylab='Intervalo de Confian?a') abline(h=12,lty=3) segments(1:nsim,li,1:nsim,ls) ############################################################################## #TEOREMA CENTRAL DO LIMITE ############################################################################## dat=c(9,1,0,7,5,6,9,5,8,8,1,0,5,7,6,5,0,2,1,2,1,8,8,8,5,2,4,8,3,1,6,5,5,7,4,1,7,3,3,3,2,8,1,8,5,8,4,0,1,9,2,1,6,9,4,4,7,6,1,7,1,9,7,9,7,2,7,7,0,8,1,6,3,8,0,5,7,4,8,6,7,0,2,8,8,7,2,5,4,1,8,6,8,3,5,8,2,7,2,4) jpeg(file='h1.jpg') hist(dat,breaks=seq(0,10,by=1),include.lowest = TRUE, right = FALSE,main='',xlab='x') dev.off() medias=NULL for (i in 1:20){ x=dat[sample(1:length(dat),size=4,replace=TRUE)] medias=c(medias,mean(x)) } jpeg(file='h2.jpg') hist(medias,breaks=seq(0,10,by=1),include.lowest = TRUE, right = FALSE,main='',xlab='m?dias') dev.off() ########## antes=rnorm(10,20,1) depois=antes+rnorm(10,2,1) d=depois-antes cbind(antes,depois,d) var(d) var(antes)+var(depois) r=cor(antes,depois) plot(antes,depois) # NOTA: var(d)=var(antes)+var(depois)-2*r*sd(antes)*sd(depois) var(antes)+var(depois)-2*r*sd(antes)*sd(depois)