En el frejol con pescado con ensalada y arroz (FPEA) se encuentran una serie de nutrientes, de los cuales en este examen se consideran los siguientes: En promedio 100 g de frejol con pescado con ensalada y arroz (FPEA) contiene 258 kcal de energia, 25.1 g de agua, 0.50 g de fibra dietaria, 23.8 g de carbohidratos totales, 4.0 g de proteinas, 16.3 g de grasa total y 0.4 g de cenizas

  1. Se analizaron 98 porciones de 100g de FPEA y en 34 de ellas el contenido de energia fue mayor de 270 kcal. Estime con un intervalo de 94% de confianza a la verdadera proporción de porciones que tienen un contenido de energia mayor de 270 keal, considerando que \(p_1 =0.02\) y \(p_2=0.04\).

\(n = 98\)

\(p = 34/98\)

\(p_1 = 0.02\)

\(p_2 = 0.04\)

\(X_i ~ Binomial(n=98,p)\)

\(T = \sum_{t=1}^{n} X_t\) es una estadistica suficiente de p

\(t_o = \sum_{t=1}^{98} X_t = 34\)

\(T = \sum_{t=1}^{98} ~ Binomial(n=98,p) \rightarrow f_T(t)= \binom{98}{t}p^t(1-p)^{20-t}I_{{0,1,..,98}}(t)\)

\(p_1 = 0.02 =\sum_{t=0}^{t_o=k_1(p)} f_T(t)=\sum_{t=0}^{34} \binom{98}{p} p^t (1-p)^{98-t} \rightarrow p =0.24964\)

\(p_2 = 0.04=\sum_{t=t_o=k_2(p)}^{98} f_T(t)=\sum_{t=0}^{34} \binom{98}{p} p^t (1-p)^{98-t} \rightarrow p = 0.43922\)

Hallando los valores por tanteo

rl1<-0.24964
rl2<-0.43922
p1_c<-1-pbinom(33,98,rl1)
p2_c<-pbinom(34,98,rl2)
c(p1_c,p2_c) #Corroborando los valores de p1 y p2
## [1] 0.02000464 0.04000189
c(rl1,rl2) #Intervalo
## [1] 0.24964 0.43922

Utilizando la generalizacion brindada en clase

y<-34
n<-98
p1 <- 0.02
p2 <- 0.04
li<-(y/(n-y+1)*qf(p1,2*y,2*(n-y+1)))/
  (1+y/(n-y+1)*qf(p1,2*y,2*(n-y+1)))
ls<-((y+1)/(n-y)*qf(1-p2,2*(y+1),2*(n-y)))/ 
  (1+(y+1)/(n-y)*qf(1-p2,2*(y+1),2*(n-y)))
print(c(li,ls),digits=8)
## [1] 0.24963596 0.43922111
  1. El número de porciones de 100 g de FPEA analizados por diatiene distribución Poisson(\(\lambda\))
  1. Con una semilla de 40 seleccione una muestra aleatoria de tamafio 28 de una distribución Poisson (10.2) y con un nivel de significacién de 0.05 concluya respecto al siguiente sistema de hipétesis UMP: \(H_o:\lambda \leq7\) o \(\lambda \geq 9\) \(H_1: 7<\lambda<9\). Presente e interprete la curva de potencia.
rm(list=ls())

# H0: p<=p1 o p>=p2
# H1: p1<p<p2
p1=7;p2=9

## Alpha=0.05, Tamano de muestra=n
a=0.05;n=5

Crear un data frame con todas las posibles combinaciones para c1<c2

cr=expand.grid('c1'=0:n,'c2'=0:n)
cr=cr[cr[,1]<cr[2],]

## debe ver la teoría
cr[,'pr1']=ppois(cr[,2]-1,p1)-ppois(cr[,1],p1)
cr[,'pr2']=dpois(cr[,1],p1)
cr[,'pr3']=dpois(cr[,2],p1)
cr[,'pr4']=ppois(cr[,2]-1,p2)-ppois(cr[,1],p2)
cr[,'pr5']=dpois(cr[,1],p2)
cr[,'pr6']=dpois(cr[,2],p2)
cr[,'g1']=a/cr[,'pr2']-cr[,'pr1']/cr[,'pr2']-(cr[,'pr3']/cr[,'pr2'])*(a*cr[,'pr2']
                                                                      +cr[,'pr1']*cr[,'pr5']-a*cr[,'pr5']-cr[,'pr2']*cr[,'pr4'])/(cr[,'pr2']*cr[,'pr6']-
                                                                                                                                    cr[,'pr3']*cr[,'pr5'])
cr[,'g2']=(a*cr[,'pr2']+cr[,'pr1']*cr[,'pr5']-a*cr[,'pr5']-cr[,'pr2']*cr[,'pr4'])/
  (cr[,'pr2']*cr[,'pr6']-cr[,'pr3']*cr[,'pr5'])
whr1=0<=cr[,'g1'] & cr[,'g1']<=1
whr2=0<=cr[,'g2'] & cr[,'g2']<=1
cr=cr[whr1 & whr2,]

comprobando si esto es igual a alpha

cr[,'a1']=cr[,'pr1']+cr[,'g1']*cr[,'pr2']+cr[,'g2']*cr[,'pr3']
cr[,'a2']=cr[,'pr4']+cr[,'g1']*cr[,'pr5']+cr[,'g2']*cr[,'pr6']
cr
##  [1] c1  c2  pr1 pr2 pr3 pr4 pr5 pr6 g1  g2  a1  a2 
## <0 rows> (or 0-length row.names)

curva de potencia

c1=cr[,'c1'];c2=cr[,'c2']
g1=cr[,'g1'];g2=cr[,'g2']

No funciono :/

  1. Con una semilla de 40 seleccione una muestra aleatoria de tamafio 28 de una distribución Poisson (10.2) y con un nivel de significación de 0.05 concluya respecto al siguiente sistema de hipétesis IUMP: \(H_0:\lambda=10.5\) vs \(H_1:\lambda \neq 10.5\). Presente e interprete la curva de potencia.

$f(x)= \(Y=\sum_{i=1}^{28}X_i\) ~ \(Poisson(28\times10.2=285.6)\)

rm(list = ls())

# H0: l=l0
# H1: l!=l0

l0=10.5

# Alpha=0.05, Sample size=n
a=0.05
n=28

#Crear un data frame con todas las posibles combinaciones para c1<c2
cr=expand.grid('c1'=0:qpois(a-0.01,n*l0),'c2'=0:qpois(1-0.01,n*l0))
cr=cr[cr[,1]<cr[,2],]

#Ver la teoría
cr[,'k1']=a-ppois(cr[,'c1']-1,n*l0)-ppois(cr[,'c2'],n*l0,lower.tail = F)
cr[,'k2']=dpois(cr[,'c1'],n*l0)
cr[,'k3']=dpois(cr[,'c2'],n*l0)
cr[,'k4']=a-ppois(cr[,'c1']-2,n*l0)-ppois(cr[,'c2']-1,n*l0,lower.tail = F)
cr[,'k5']=dpois(cr[,'c1']-1,n*l0)
cr[,'k6']=dpois(cr[,'c2']-1,n*l0)

cr[,'g2']=(cr[,'k4']-cr[,'k1']*cr[,'k5']/cr[,'k2'])/(cr[,'k6']-cr[,'k3']*cr[,'k5']/cr[,'k2'])
cr[,'g1']=(cr[,'k1']-cr[,'g2']*cr[,'k3'])/cr[,'k2']

whr1 = 0<=cr[,'g1'] & cr[,'g1']<=1
whr2 = 0<=cr[,'g2'] & cr[,'g2']<=1
cr=cr[whr1 & whr2,]

#chequeando para ver si alpha es igual
cr[,'a1']=ppois(cr[,'c1']-1,n*l0)+ppois(cr[,'c2'],n*l0,lower.tail = F)+
  cr[,'g1']*dpois(cr[,'c1'],n*l0)+cr[,'g2']*dpois(cr[,'c2'],n*l0)
cr[,c('c1','c2','g1','g2','a1')]
##        c1  c2       g1        g2   a1
## 87182 261 328 0.521422 0.2434018 0.05
#curva de potencia
pwr=function(x) ppois(cr[,'c1']-1,n*x)+ppois(cr[,'c2'],n*x,lower.tail = F)+
  cr[,'g1']*dpois(cr[,'c1'],n*x)+cr[,'g2']*dpois(cr[,'c2'],n*x)
curve(pwr,main='Potencia para prueba IUMP de dos colas: Poisson',xlab='Valores de lambda'
      ,ylab='Probabilidad de rechazar H0',xlim=c(0,4),ylim=c(0,1),
      col='green2',lwd=3)
abline(h=a,v=l0,lty=2)

#simulacion co la computadora
c1=cr[,'c1']
c2=cr[,'c2']
g1=cr[,'g1']
g2=cr[,'g2']
loops=1000000
x=rpois(loops,n*10.5)
tmp=(x<c1 | c2<x) | ((runif(loops)<=g1)&(x==c1)) | ((runif(loops)<=g2)&(x==c2))
(tmp=table(tmp))
## tmp
##  FALSE   TRUE 
## 950333  49667
tmp[2]/loops
##     TRUE 
## 0.049667
segments(x0=10.5,y0=0,y1=tmp[2]/loops,col='red')
segments(x0=0,y0=tmp[2]/loops,x1=10.5,col='red')

# "colas iguales" enfoque insesgado aproximado

cr[,'c1b']=qpois(a/2,n*l0)
cr[,'c2b']=qpois(a/2,n*l0,lower.tail = F)
cr[,'g1b']=(a/2-ppois(cr[,'c1b']-1,n*l0))/dpois(cr[,'c1b'],n*l0)
cr[,'g2b']=(a/2-ppois(cr[,'c2b'],n*l0,lower.tail = F))/dpois(cr[,'c2b'],n*l0)

cr[,'a1b']=ppois(cr[,'c1b']-1,n*l0)+ppois(cr[,'c2b'],n*l0,lower.tail=F)+
  cr[,'g1b']*dpois(cr[,'c1b'],n*l0)+cr[,'g2b']*dpois(cr[,'c2b'],n*l0)
cr[,c('c1','c2','g1','g2','a1','g1b','g2b','a1b')]
##        c1  c2       g1        g2   a1       g1b       g2b  a1b
## 87182 261 328 0.521422 0.2434018 0.05 0.3630022 0.4155272 0.05
equal.tails=function(x) ppois(cr[,'c1b']-1,n*x)+ppois(cr[,'c2b'],n*x,lower.tail = F)+
  cr[,'g1b']*dpois(cr[,'c1b'],n*x)+cr[,'g2b']*dpois(cr[,'c2b'],n*x)
  1. El contenido de proteinas en 100 g de FPEA tiene distribución Normal(\(\mu =4\),\(\sigma^2=1.2\)) Con una semilla de 42 seleccione una muestra aleatoria detamafio 40 de una distribución Normal (\(\mu=4\),\(\sigma=1.2\)) y con un nivel de significación de 0.05 concluya respecto al siguiente sistema de hipotesis IUMP: \(H_0: =1.19\) vs \(H_1: \sigma^2 \neq 1.19\). Calcule la potencia de la prucba cuando en realidad \(\sigma^2 =1.52\)

\(f(x)=\frac{1}{\sqrt{2\pi}\sigma}[exp(-\frac{x^2}{2\sigma^2})]=\frac{exp[tx^2]}{\sqrt{2\pi}\sigma}\)

Aplicando la función de verosimilitud, la productoria. Del cual podemos decir que:

\(Y=\sum_{t=1}^{n}x^2_t\) es suficiente para \(\sigma^2\)

\(\frac{Y}{\sigma^2} = \chi^2(n)\). \(Y\) ~ \(\sigma^2 \chi^2(n)\)

Notamos que \(-\frac{1}{2\sigma^2}\) es una función 1-1 de \(\sigma^2 > 0\)

\(\phi(\underline{x})=\left\{\begin{matrix} 1 & si\ Y <c_1 \ o \ c_2<Y \ (c_1,c_2) \\ 0 & si\ c_1 <Y<c_2 \end{matrix}\right.\)

\(W=\frac{Y}{\sigma^2_0}\)

Para una hipotesis simple vs compuesta se cumple lo siguiente:

(1)…

\[E_{\sigma^2_0}[\phi(x)]=\alpha\]

(2)…

\[E_{\sigma^2_0}[W\phi(\underline{x})]=\alpha E_{\sigma^2_0}[W]\] Ello nos quedarian dos ecuaciones con dos incognitas, por lo cual preferimos hacer una aproximacion de colas iguales.

\(\alpha = \alpha E_{\sigma^2_0}[\phi(\underline{x})]=P_{\sigma^2_0}(W<c_1)+P_{\sigma^2_0}(c_2'<W)\)

Se establece: \(P_{\sigma^2_0}(W<c_1)=P_{\sigma^2_0}(c_2'(c_2<W) = \frac{\alpha}{2}\)

\(c_1' = \chi^2_{\frac{\alpha}{2}}(n) \to c_1 =\sigma ^2_0\chi^2_{\frac{\alpha}{2}}(n)\)

\(c_2' = \chi^2_{1-\frac{\alpha}{2}}(n) \to c_2 =\sigma ^2_0\chi^2_{1-\frac{\alpha}{2}}(n)\)

rm(list = ls())

#Planteo de la hipotesis:

# H0: Sigma^2=Sigma0^2
# H1: Sigma^2!=Sigma0^2

Datos

s2=1.19  ## Sigma0^2 

# Alpha=0.05, Sample size=n
n=40;a=0.05

Elaboramos un sistema de ecuaciones

alpha.left.tail=seq(0.001,a-0.001,0.001)
c1=qchisq(alpha.left.tail,n)
c2=qchisq(1-(a-alpha.left.tail),n)

c1b=qchisq(alpha.left.tail,n+2)
c2b=qchisq(1-(a-alpha.left.tail),n+2)

Gráfico para ver lo que se esta resolviendo.

\(c_1\) y \(c_2\) es la linea negra \(c_1'\) y \(c_2'\) es la linea roja

plot(c1,c2,,type='l')
lines(c1b,c2b,col='red')

cr=function(c1) qchisq(1-(a-pchisq(c1,n)),n)-qchisq(1-(a-pchisq(c1,n+2)),n+2)

## use uniroot to find c1
crit.reg1=uniroot(cr,range(c1),tol=10^(-15))$root

## usar c1 para hallar c2
crit.reg2=qchisq(1-(a-pchisq(crit.reg1,n)),n)
abline(v=crit.reg1,h=crit.reg2)

Región critica

c("c1'"=crit.reg1,"c2'"=crit.reg2,'c1'=s2*crit.reg1,'c2'=s2*crit.reg2)
##      c1'      c2'       c1       c2 
## 24.87866 60.27483 29.60560 71.72705

Teniendo como resultado:

Funcion potencia

La probabilidad de rechazar la \(H_o\)

pwr=function(s) pchisq(s2*crit.reg1/s,n)+pchisq(s2*crit.reg2/s,n,lower.tail=F)
curve(pwr,xlim=c(0,2*s2),main='Potencia para una prueba IUMP de 2 colas: Varianza Normal',
      xlab='Varianza Verdadera',ylab='Probabilidad de rechazar H0',ylim=c(0,1),
      col='green2',lwd=3)
abline(h=a,v=s2,lty=2,col='purple')

#computer simulation
loops=500000
sig2=1.52

set.seed(42)
x=matrix(rnorm(n*loops,mean=0,sd=sqrt(sig2)),nrow=loops,ncol=n)
x=rowSums(x^2)
tmp=x<s2*crit.reg1 | s2*crit.reg2<x
(tmp=table(tmp))
## tmp
##  FALSE   TRUE 
## 396942 103058
tmp[2]/loops
##     TRUE 
## 0.206116
segments(x0=sig2,y0=0,y1=tmp[2]/loops,col='red')
segments(x0=0,y0=tmp[2]/loops,x1=sig2,col='red')

Podemos observar que medida que se aleja de 1.19 aumenta la potencia de la prueba.

Una prueba IUMP es aquella cuya esperanza es menor igual que alfa cuando estamos en la hipotesis nula y esa esperanza es mayor igual que alfa cuando estamos en la hipotesis alterna.

En donde el error tipo I es como maximo el valor de alfa, es decir, 0.05 porque es el valor minimo de la potencia.

Finalmente para un valor de la varianza igual a 1.52 la probabilidad de rechazar la hipotesis nula (\(H_o\)) es igual a 0.206