Universidade Federal do Paraná
Curso de Estatística
CE 083 -
Estatística Computacional I - 2014/2
Prof. Dr. Walmes Marques Zeviani
##-----------------------------------------------------------------------------
## Só para curiosos. Ferramentas interativas com Shiny.
browseURL("http://shiny.albany.edu/stat/binomial/")
browseURL("https://econometricsbysimulation.shinyapps.io/OLS-App/")
##-----------------------------------------------------------------------------
## Definições da sessão.
## Pacote para janelas interativas.
require(rpanel)
##-----------------------------------------------------------------------------
## Distribuição binomial.
pb <- function(panel){
with(panel,
{
x <- 0:size
px <- dbinom(x, size=size, prob=prob)
plot(px~x, type="h", ylim=c(0,max(c(px), 0.5)))
if(showEX){
abline(v=size*prob, col=2)
}
})
panel
}
## args(rp.slider)
## help(rp.slider, help_type="html")
## help(rp.control, help_type="html")
panel <- rp.control(title="Binomial", size=c(300,100))
rp.slider(panel, size, from=2, to=80, initval=10, resolution=1,
action=pb, showvalue=TRUE, title="size")
rp.slider(panel, prob, from=0.01, to=0.99, initval=0.5, resolution=0.01,
action=pb, showvalue=TRUE, title="prob")
rp.checkbox(panel, showEX, action=pb, title="E(X)",
labels="Mostrar o valor esperado?")
##-----------------------------------------------------------------------------
## Distribuição de Poisson.
pp <- function(panel){
with(panel,
{
x <- 0:100
px <- dpois(x, lambda=lambda)
plot(px~x, type="h", ylim=c(0,max(c(px),0.5)))
if(showEX){
abline(v=lambda, col=2)
}
})
panel
}
panel <- rp.control(title="Poisson", size=c(300,100))
rp.slider(panel, lambda, from=0.5, to=90, initval=10, resolution=0.25,
action=pp, showvalue=TRUE, title="lambda")
rp.checkbox(panel, showEX, action=pp, title="E(X)",
labels="Mostrar o valor esperado?")
##-----------------------------------------------------------------------------
## Distribuição binomial negativa.
pbn <- function(panel){
with(panel,
{
x <- 0:size
px <- dnbinom(x, size=size, prob=prob)
plot(px~x, type="h", ylim=c(0,max(c(px), 0.5)))
if(showEX){
abline(v=size*(1-prob)/prob, col=2)
}
})
panel
}
panel <- rp.control(title="Binomial negativa", size=c(300,100))
rp.slider(panel, size, from=2, to=80, initval=10, resolution=1,
action=pbn, showvalue=TRUE, title="size")
rp.slider(panel, prob, from=0.01, to=0.99, initval=0.5, resolution=0.01,
action=pbn, showvalue=TRUE, title="prob")
rp.checkbox(panel, showEX, action=pbn, title="E(X)",
labels="Mostrar o valor esperado?")
##-----------------------------------------------------------------------------
## Distribuição hipergeométrica.
## m the number of white balls in the urn.
## n the number of black balls in the urn.
## k the number of balls drawn from the urn.
## x the number of white balls drawn without replacement.
ph <- function(panel){
with(panel,
{
x <- max(c(0, k-m)):min(c(k, m))
## p(x) = choose(m, x) choose(n, k-x) / choose(m+n, k)
px <- dhyper(x, m=m, n=n, k=k)
plot(px~x, type="h", ylim=c(0, max(c(px), 0.5)))
if(showEX){
abline(v=k*m/(m+n), col=2)
}
})
panel
}
panel <- rp.control(title="Hipergeométrica", size=c(300,100))
rp.slider(panel, m, from=5, to=30, initval=10, resolution=1,
action=ph, showvalue=TRUE, title="Brancas")
rp.slider(panel, n, from=2, to=15, initval=5, resolution=1,
action=ph, showvalue=TRUE, title="Pretas")
rp.slider(panel, k, from=2, to=15, initval=5, resolution=1,
action=ph, showvalue=TRUE, title="Retiradas")
rp.checkbox(panel, showEX, action=ph, title="E(X)",
labels="Mostrar o valor esperado?")
##-----------------------------------------------------------------------------
## Distribuição normal.
pn <- function(panel){
with(panel,
{
curve(dnorm(x, mean=mean, sd=sd), -5, 5,
ylim=c(0,1))
if(showEX){
abline(v=mean, col=2)
}
if(showVX){
d <- dnorm(mean+sd, mean, sd)
segments(mean-sd, d, mean+sd, d, col=2)
}
})
panel
}
panel <- rp.control(title="Normal", size=c(300,100))
rp.slider(panel, mean, from=-4, to=4, initval=0, resolution=0.1,
action=pn, showvalue=TRUE, title="mean")
rp.slider(panel, sd, from=0, to=3, initval=1, resolution=0.1,
action=pn, showvalue=TRUE, title="sd")
rp.checkbox(panel, showEX, action=pn, title="E(X)",
labels="Mostrar o valor esperado?")
rp.checkbox(panel, showVX, action=pn, title="sd(X)",
labels="Mostrar o desvio-padrão?")
##-----------------------------------------------------------------------------
## Distribuição beta.
pbt <- function(panel){
with(panel,
{
curve(dbeta(x, shape1=exp(sh1), shape2=exp(sh2)), 0, 1,
ylim=c(0,7))
if(showEX){
abline(v=exp(sh1)/(exp(sh1)+exp(sh2)), col=2)
}
})
panel
}
panel <- rp.control(title="Beta", size=c(300,100))
rp.slider(panel, sh1, from=-5, to=5, initval=0, resolution=0.1,
action=pbt, showvalue=TRUE, title="log(shape1)")
rp.slider(panel, sh2, from=-5, to=5, initval=0, resolution=0.1,
action=pbt, showvalue=TRUE, title="log(shape2)")
rp.checkbox(panel, showEX, action=pbt, title="E(X)",
labels="Mostrar o valor esperado?")
##-----------------------------------------------------------------------------
## Distribuição gamma.
pg <- function(panel){
with(panel,
{
curve(dgamma(x, shape=shape, scale=scale), 0, 50,
ylim=c(0, 0.25))
if(showEX){
abline(v=shape*scale, col=2)
}
})
panel
}
panel <- rp.control(title="Gamma", size=c(300,100))
rp.slider(panel, shape, from=0.1, to=20, initval=5, resolution=0.1,
action=pg, showvalue=TRUE, title="shape")
rp.slider(panel, scale, from=0.1, to=10, initval=3, resolution=0.1,
action=pg, showvalue=TRUE, title="scale")
rp.checkbox(panel, showEX, action=pg, title="E(X)",
labels="Mostrar o valor esperado?")
##-----------------------------------------------------------------------------
## Distribuição Weibull.
pw <- function(panel){
with(panel,
{
curve(dweibull(x, shape=shape, scale=scale), 0, 50,
ylim=c(0, 0.25))
if(showEX){
abline(v=scale*gamma(1+1/shape), col=2)
}
})
panel
}
panel <- rp.control(title="Weibull", size=c(300,100))
rp.slider(panel, shape, from=0.1, to=10, initval=5, resolution=0.1,
action=pw, showvalue=TRUE, title="shape")
rp.slider(panel, scale, from=0.1, to=30, initval=20, resolution=0.1,
action=pw, showvalue=TRUE, title="scale")
rp.checkbox(panel, showEX, action=pw, title="E(X)",
labels="Mostrar o valor esperado?")
##-----------------------------------------------------------------------------
## Função que calcula probabilidades para uma v.a. X ~ Binomial(n, p).
pbin <- function(x, n, p){
choose(n,x)*p^x*(1-p)^(n-x)
}
choo
## Error: object 'choo' not found
## Pr(X=20 | n=25, p=0.9): probabilidade em um ponto.
pbin(20, n=25, p=0.9)
## Pr(20<=X<=25 | n=25, p=0.9): probabilidade em um intervalo.
sum(pbin(20:25, n=25, p=0.9))
## Funções do R para fazer o mesmo.
dbinom(20, size=25, prob=0.9) ## Pr(X=20).
pbinom(20, size=25, prob=0.9) ## Pr(X<=20).
1-pbinom(19, size=25, prob=0.9) ## 1-Pr(X<=19) = Pr(X>19) = Pr(X>=20).
pbinom(19, size=25, prob=0.9,
lower.tail=FALSE) ## Pr(X>19).
##-----------------------------------------------------------------------------
## Para a distribuição normal.
## Gera números aleatórios.
rnorm(10, mean=1.8, sd=0.1)
## Gráfico da densidade.
curve((1/sqrt(2*pi*0.1^2)*exp(-(x-1.8)^2/(2*0.1^2))), 1.3, 2.3)
curve(dnorm(x, mean=1.8, sd=0.1), add=TRUE, col=2, lty=2, lwd=2)
## Gráfico da densidade integrada ou probabilidade acumulada.
curve(pnorm(x, mean=1.8, sd=0.1), 1.3, 2.3)
## Pr(1.75<X<=1.95).
pnorm(1.95, 1.8, 0.1)-pnorm(1.75, 1.8, 0.1)
## Representação da probabilidade correspondente ao evento.
xseq <- seq(1.75, 1.95, length.out=50)
fseq <- dnorm(xseq, mean=1.8, sd=0.1)
curve(dnorm(x, mean=1.8, sd=0.1), 1.3, 2.3)
polygon(x=c(head(xseq, 1), xseq, tail(xseq, 1)),
y=c(0, fseq, 0), col="gray90")
curve(pnorm(x, mean=1.8, sd=0.1), 1.3, 2.3)
abline(v=c(1.75, 1.95), h=pnorm(c(1.75, 1.95), mean=1.8, sd=0.1), lty=2)
##-----------------------------------------------------------------------------
## A função quantil em distribuições discretas.
## Definição:
## quantil(p) = inf{x \in R : p <= F(x)}
q <- qpois(0.4, lambda=10); q
curve(ppois(x, lambda=10), 0, 25, type="s")
abline(h=0.40, v=q, lty=2)
