O conceito de funcao densidade de probabilidade se assemelha ao conceito de funcao de probabilidade, que serve para o caso de variaveis aleatorias discretas.
Neste caso uma variavel aleatoria discreta tem um numero de possiveis ocorrencias. Aqui, a variavel aleatoria resultado do lancamento da moeda tem apenas 2 possiveis ocorrencias: cara ou coroa. Por isso, a funcao de probabilidade a ela associada tambem so podera assumir 2 valores, 0.4 e 0.6 (que neste caso a moeda esta viciada). Assim:
x <- 1:4
fx <- dbinom(x, 4, 0.4)
plot(x, fx, type ='h', main = "Funcao de distribuicao de probabilidade", lwd = 5, col = "blue4")
Fx <- pbinom(x, 4, 0.4)
plot(x, Fx, type='S', main = "Funcao de distribuicao acumulada", lwd = 5, col = "magenta")
neste caso sendo a binomial uma distribuicao discreta de probabilidades. Os graficos obtidos sao diferentes da distribuicao continuas, como a distribuicao normal. A funcao de probabilidade seria uma descricao das probabilidades associadas com os valores possiveis de X, enquanto a funcao distribuicao acumulada calcula a probabilidade acumulada para o valor de X. Para uma variavel discreta, a distribuicao seria apenas uma lista de valores possiveis, com as suas probabilidades associadas
fx
## [1] 0.3456 0.3456 0.1536 0.0256
fx_h <- dbinom(x, 4, 0.5)
fx_h
## [1] 0.2500 0.3750 0.2500 0.0625
Assim, percebe-se que as probabilidades sao diferentes. Isto se da porque a moeda esta viciada.
Do exercicio podemos resumir os dados observados em numeros de frequencias da incidencia da doenca. De costume consideramos uma tabela de distribuicao de frequencias com \(k\) categorias tal que essas frequencias devam ser superiores ou iguais a \(2\). Assim:
incid <- c(0,1,2,3,4,5)
observadas <- c(4,4,6,6,2,4)
esperadas <- c(2.6,2.6,7.8,7.8,2.6,2.6)
A tabela gerada temos que:
tabela <- data.frame(incid,observadas,esperadas)
fator_incid <- as.factor(incid)
tabela$incid <- fator_incid
tabela
## incid observadas esperadas
## 1 0 4 2.6
## 2 1 4 2.6
## 3 2 6 7.8
## 4 3 6 7.8
## 5 4 2 2.6
## 6 5 4 2.6
Existem testes de adequabilidade para saber se a amostra que temos segue um possivel modelo probabilistico. Estes testes de adequabilidade sao chamados de testes de aderencia.
Considere a probabilidade especificada \(p_i\) a categoria \(i\) de acordo com o modelo que desejamos verificar, \(i=1,.,k\). Neste caso ha \(6\) categorias que desejamos testar (numero de casos de cancer em parentes proximos). As hipoteses deste teste sao:
\[H_0 : p_0 = 1/10, p_1 = 1/10,..., p_5 = 1/10 \] \[H_1 : Existe\ pelo\ menos\ uma\ diferente.\]
chisq.test(tabela$observadas,
p = c(1/10, 1/10, 3/10, 3/10, 1/10, 1/10),
correct = FALSE)
## Warning in chisq.test(tabela$observadas, p = c(1/10, 1/10, 3/10, 3/10, 1/10, :
## Chi-squared approximation may be incorrect
##
## Chi-squared test for given probabilities
##
## data: tabela$observadas
## X-squared = 3.2308, df = 5, p-value = 0.6645
A funcao chisq.test do teste mostra que, com base no nivel de significancia de \(5\%\) nao rejeita-se \(H_0\) pois o \(p-valor\) do teste seria maior que o nivel de significancia \(0.05\). Ou seja, nao foi encontrado evidencias nos dados para desacreditar do modelo proposto.
O grafico pode ser feito por
plot(tabela$observadas, main = "frequencias observadas")
plot(tabela$esperadas, main = "frequencias teoricas")
Os dados observados parecem seguir a distribuicao proposta teorica, no entanto ha uma discrepancia nestes valores. A partir do grafico nao parece haver uma concordancia entre a amostra e o modelo teorico. Alem disso o tamanho da amostra influencia no resultado.
Sera calculado nos, exercicios 3) e 4), as probabilidades associadas a cada evento do espaco amostral um valor para a variavel aleatoria X ou suas probabilidades acumuladas. Assim:
1-pbinom(13,15,0.4)
## [1] 2.523293e-05
pbinom(10, 15, 0.4) - pbinom(8, 15, 0.4)
## [1] 0.08569975
pbinom(1, 15, 0.4) + (1 - pbinom(10, 15, 0.4))
## [1] 0.0145197
(1-pbinom(10, 15, 0.4)) + (1-pbinom(13, 15, 0.4))
## [1] 0.009372894
c <- pbinom(6, 15, 0.4)
d <- 1-pbinom(3, 15, 0.4)
c*d
## [1] 0.5546239
a<- pbinom(14, 15, 0.4) - pbinom(10, 15, 0.4)
b<- 1-pbinom(10, 15, 0.4)
a/b
## [1] 0.9998851
pnorm(116, mean = 90, sd = 100)
## [1] 0.6025681
pnorm(79, mean = 90, sd = 100)
## [1] 0.4562047
pnorm(76, mean = 90, sd = 100)
## [1] 0.44433
pnorm(111, mean = 90, sd = 100) - pnorm(84, mean = 90, sd = 100)
## [1] 0.1070883
pnorm(101, mean = 90, sd = 100) - pnorm(79, mean = 90, sd = 100)
## [1] 0.08759063
A funcao de densidade descreve a verossimilhanca de uma variavel aleatoria tomar um valor dado. Ou seja, a probabilidade da variavel aleatoria cair em uma faixa particular esta sob a regiao da area do grafico.
#-----------
# lambda: 5
#-----------
lambda <- 5
plot(dpois(x, lambda), type = "h", lwd = 2,
main = "funcao de distribuicao de densidade",
ylab = "P(X = x)", xlab = "numero de eventos")
# Legenda
legend("topright", legend = c("5"),
title = expression(lambda), title.adj = 0.75,
lty = 1, col = 1:3, lwd = 2, box.lty = 0)
plot(function(x) dnorm(x, 90, 100), 60, 140, ylab='f(x)')
title('Distribuicao Normal X ~ N(80, 100)')
#tamanho do eixo x
N=10000
y_rnorm <- rnorm(N, mean = 90, sd = 100)
y_rnorm2 <- rnorm(N, mean = 90, sd = 80)
y_rnorm3 <- rnorm(N, mean = 85, sd = 100)
plot(density(y_rnorm),
xlim = c(- 400, 600),
ylim = c(0, 0.0055),
main = "densidade da distribuicao normal", lwd = 2)
lines(density(y_rnorm2), col = "coral2", lwd = 2)
lines(density(y_rnorm3), col = "green3", lwd = 2)
legend("topright",
legend = c("Media = 100; SD = 64",
"Media = 90; SD = 80",
"Media = 85; SD = 100"),
col = c("black", "coral2", "green3"),
lty = 1)
library(stats)
#tamanho do eixo x
z <- 0:100
#-----------
# gl: 1
#-----------
gl <- 1
plot(dchisq(z, gl,ncp = 0, log = FALSE), type = "h", lwd = 2,
main = "funcao de distribuicao de densidade",
ylab = "P(X = x)", xlab = "numero de eventos")
#-----------
# gl: 2
#-----------
gl <- 2
lines(dchisq(z, gl,ncp = 0, log = FALSE), type = "h", lwd = 2, col = rgb(1,0,0, 0.7))
#-----------
# gl: 5
#-----------
gl <- 5
lines(dchisq(z, gl,ncp = 0, log = FALSE), type = "h", lwd = 2, col = rgb(0, 1, 0, 0.7))
# legenda
legend("topright", legend = c("1", "2", "5"),
title = "graus de liberdade (gl)", title.adj = 0.75,
lty = 1, col = 1:3, lwd = 3, box.lty = 1)