##### Skript pro 8. cvičení ##### ----------------------------------- ## ## téma: binomické rozdělení, normální rozdělení ## interval spolehlivosti pro střední hodnotu N(.,.) ## getwd() # měli bychom dostat "J:/statistika", pokud ne, pak: setwd("J:/statistika") getwd() # měli bychom dostat "J:/statistika" rm(list=ls()) ## ## # Ke zkoušce se daný den dostavilo 25 studentů, o kterých předpokládáme, # že jsou stejně dobře připraveni, pravděpodobnost neúspěchu je 15 %. # Které rozdělení má náhodná veličina Y (počet neúspěšných studentů)? # Který počet neúspěšných je nejpravděpodobnější? cbind("k"=0:10,"P(Y=k)"=round(dbinom(0:10,25,0.15),6)) cbind("k"=k<-0:10,"P(Y=k)"=round(dbinom(k,25,0.15),6)) # Který z jevů # A = všichni zkoušku udělají, # B - neprojde právě sedm studentů # je méně pravděpodobný? # dbinom(c(0,7),25,0.15) ## library(TeachingDemos) vis.binom() ## # Odpovězte na podobné otátky, ale pro případ 75 studentů, # A - neprojde pět studentů # B - neprojde šestnáct stdudentů. # Porovnejte přesný výpočet s normální aproximací. # dbinom(c(5,16),75,0.15) dnorm(c(5,16),75*0.15,sqrt(75*0.15*0.85)) cbind(0:20,dbinom(0:20,75,0.15),dnorm(0:20,75*0.15,sqrt(75*0.15*0.85))) cbind("k"=k<-0:20,"binomické"=round(dbinom(k,75,0.15),6), "norm. aprox."=round(dnorm(k,75*0.15,sqrt(75*0.15*0.85)),6)) ## vis.normal() ## # Předpokládejme, že výšky dospělých mužů produktivního věku v cm lze popsat # normálním rozdělením N(180, 50). Jaká část této populace má výšku # v rozmezí od 170 do 190 cm? pnorm(190,180,sqrt(50))-pnorm(170,180,sqrt(50)) # U jaké části populace takovou výšku NAMĚŘÍME ? pnorm(190.5,180,sqrt(50))-pnorm(169.5,180,sqrt(50)) ## ## normální rozdělení vzniká jako součet velkého počtu ## nezávislých stejně rozdělených NV ## data(Matky) ls() length(matky) summary(matky) hist(matky) (mu = mean(matky)) (sigma2 = var(matky)) table(matky) sample(matky,15) sample(matky,15) sample(matky,15) mean(sample(matky,15)) mean(sample(matky,15)) mean(sample(matky,15)) # B = 100 n = 10 průměry = rep(NA,B) for (b in 1:B) průměry[b] = mean(sample(matky,n)) # c("odhad"=mean(průměry),"skuečnost"=mu) c("odhad"=var(průměry),"teorie"=sigma2/n) hist(průměry) # vyzkoušet i pro jiná n # # interval spolehlivosti (konfidenční interval) n = 15 (x = sample(matky,n)) mean(x) sd(x) (n = length(x)) sd(x)/sqrt(n) qt(0.975,n-1)*sd(x)/sqrt(n) (tKrit = qt(0.975,n-1)) c(-1,1)*tKrit*sd(x)/sqrt(n) mean(x)+c(-1,1)*tKrit*sd(x)/sqrt(n) t.test(x) # # opakovaný výpočet CI pro opakované výběry průměry = rep(NA,B) CI = matrix(NA,nrow=B,ncol=2) # for (b in 1:B) { x = sample(matky,n) průměry[b] = mean(x) CI[b,] = mean(x)+c(-1,1)*tKrit*sd(x)/sqrt(n) } c("odhad"=mean(průměry),"skuečnost"=mu) c("odhad"=var(průměry),"teorie"=sigma2/n) chyba = (CI[,1]>=mu)|(CI[,2]<=mu) data.frame(CI[,1],CI[,2],chyba) ("rel. četnost chyby" = mean(chyba)) hist(průměry) rozsah = range(CI) plot(1:B,průměry,pch=3,ylim=rozsah,xlab="výběr",ylab="věk") abline(h=mu,col="red") for (b in 1:B) arrows(b,CI[b,1],b,CI[b,2],length=0,col=ifelse(chyba[b],"red","blue")) ## ######################################################################## ## data(Kojeni) names(Kojeni) with(Kojeni,t.test(vyskaM)) with(Kojeni,t.test(vyskaO)) ## Co říkají právě spočítané konfidenční intervaly? with(Kojeni,t.test(vyskaM),conf.level = 0.99) with(Kojeni,t.test(vyskaO),conf.level = 0.99) with(Kojeni,t.test(vyskaM),conf.level = 0.90) with(Kojeni,t.test(vyskaO),conf.level = 0.90) ## data(Policie) names(Policie) t.test(Policie$height) t.test(Policie$weight) ## zjistit výšky přítomných x = c( ) t.test(x) ## ## interval spolehlivosti pro podíl ## ## Plan - zda těhotenství označeno za plánované # převést na 0 - 1 veličinu plan = as.numeric(Kojeni$Plan==1) summary(plan) (n = length(plan)) (Y = sum(plan)) (f = mean(plan)) # relativní četnost jedniček f+c(-1,1)*sqrt(f*(1-f)/n)*1.96 # přibližný interval spolehlivosti prop.test(Y,n) binom.test(Y,n)