R Markdown

This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.

When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:

###Tamaño de muestra para determinar proporción:

##En la fórmula del libro Ni es el tamaño de cada estrato, Pi es el valor de referencia (historial de defectuosos), qi es los que no estan defectuosos, ai es Ni sobre total de items.

##Crear función para el cálculo (estimar tamaño de muestra):


f_sam_prop <- function(Ni, pi, ti, e, z=1.96){
  N = sum(Ni)  ##tamaño población
  L= length(Ni) ##número de estratos
  qi = 1 -pi    ##Proporción no defectuosos
  ai = Ni/sum(Ni) ##peso de cada estrato
  D = (e/z)**2
  
##ti es el tiempo de operación en cada turno (en minutos)
  
  numer = sum(Ni**2 * pi*qi/ai)
  denom = (N**2*D) + sum(Ni*pi*qi)
  n= ceiling(numer/denom)
  
  
  ###Asignación proporcional por estrato
  ni = round(ai*n)  ##tamaño de muestra por estrato
  ki= floor(ti/ni)  ###tiempo de muestreo sistemático por estrato
  return(list(ni,ki))
}


f_sam_prop(Ni=c(6000, 5000, 4000),
           pi = c(0.04, 0.05, 0.08),
           ti = c(480, 480, 480),
           e= 0.05)
## [[1]]
## [1] 31 26 21
## 
## [[2]]
## [1] 15 18 22
##Asumamos que al hacer los recorridos yo obtuve un nuevo valor de proporción de defectuosos:
##Ejemplo: 

a= 2/5*0.03 + 1/3*0.04 + 4/15*0.06
(2/5*0.03)/a
## [1] 0.2903226
(1/3*0.04)/a
## [1] 0.3225806
(4/15*0.06)/a
## [1] 0.3870968
###TAREA HACER EL EJERCICIO CON SOLO 2 TURNOS
###Muestreo espacial con cobertura:
palmas = expand.grid(x=seq(0, 112, 7),
                     y=seq(0, 144, 9))
nrow(palmas)
## [1] 289
set.seed(123)
###Información auxiliar (producción previa)
p_racimo_u = rnorm(289, 17, 1.8)
p_racimo_p = rnorm(289, 17, 1.8)
CaMg_h17 = runif(289, 1.8, 2.0)  ##relación calcio magnesio hoja 17
CaMg_s = runif(289, 1.2, 1.4) ##relación calcio magnesio suelo

hibrid = rep(c("h1", "h2"), c(144, 145))



df= data.frame(palmas, 
               p_racimo_u, p_racimo_p,
               CaMg_h17, CaMg_s)
head(df)
##    x y p_racimo_u p_racimo_p CaMg_h17   CaMg_s
## 1  0 0   15.99114   17.30252 1.813162 1.350780
## 2  7 0   16.58568   19.10309 1.970874 1.212619
## 3 14 0   19.80567   18.89753 1.902082 1.288496
## 4 21 0   17.12692   19.06147 1.842000 1.300134
## 5 28 0   17.23272   15.96056 1.875169 1.254049
## 6 35 0   20.08712   20.60447 1.841473 1.319760
library(clhs)
## Warning: package 'clhs' was built under R version 4.1.3
 res <- clhs(df, size=15, progress=FALSE, simple=TRUE)
 df[res,"muestreo"] = "si"
 df$muestreo[is.na(df$muestreo)] = "no"
 
 library(ggplot2)
## Warning: package 'ggplot2' was built under R version 4.1.3
 ggplot(df)+
   aes(x,y, fill=muestreo)+
   geom_tile(color="white")

##Demostrar que el modelo se acerca bastante bien a los "datos reales"  
mean(df$p_racimo_u)
## [1] 17.02191
mean(df$p_racimo_p)
## [1] 17.09966
mean(df$p_racimo_u[df$muestreo == "si"])
## [1] 16.68358
mean(df$p_racimo_p[df$muestreo == "si"])
## [1] 16.64809
###En la funcion "size 15" no alcanzó una cobertura de todo el lote. Probar con 20 palmas:
set.seed(124)

res2 <- clhs(df, size=20, progress=FALSE, simple=TRUE)
## Warning: NAs introducidos por coerción
 df[res2,"muestreo"] = "si"
 df$muestreo[is.na(df$muestreo)] = "no"
 
 library(ggplot2)
 ggplot(df)+
   aes(x,y, fill=muestreo)+
   geom_tile(color="white")

###ENTREGA LAURA LAVERDE - MÉTODOS MULTIVARIADOS

#EJERCICIO 1
###Tamaño de muestra para determinar proporción (con solo dos turnos: diurno y nocturno)
##Crear función para el cálculo (estimar tamaño de muestra):


f_sam_prop <- function(Ni, pi, ti, e, z=1.96){
  N = sum(Ni)  ##tamaño población
  L= length(Ni) ##número de estratos
  qi = 1 -pi    ##Proporción no defectuosos
  ai = Ni/sum(Ni) ##peso de cada estrato
  D = (e/z)**2
  
##ti es el tiempo de operación en cada turno (en minutos)
  
numer = sum(Ni**2 * pi*qi/ai)
denom = (N**2*D) + sum(Ni*pi*qi)
n= ceiling(numer/denom)
  
  
###Asignación proporcional por estrato
ni = round(ai*n)  ##tamaño de muestra por estrato
ki= floor(ti/ni)  ###tiempo de muestreo sistemático por estrato
return(list(ni,ki))
}


f_sam_prop(Ni=c(8500, 6500),
           pi = c(0.04, 0.08),
           ti = c(720, 720),
           e= 0.05)
## [[1]]
## [1] 46 36
## 
## [[2]]
## [1] 15 20
###Ahora con los datos del muestreo:
###Cálculo de ai (peso de cada turno)
8500/15000 ##ai_diurno 
## [1] 0.5666667
6500/15000 ##ai_nocturno 
## [1] 0.4333333
####Si yo quiero saber cual es la PROBABILIDAD CONDICIONAL:
###ASUMAMOS que la probabilidad de defectuosos basada en muestreo en turno diurno es 0.03 y en turno nocturno es 0.06:
a= 8500/15000*0.03 + 6500/15000*0.06
(8500/15000*0.03)/a
## [1] 0.3953488
(6500/15000*0.06)/a
## [1] 0.6046512
###Concluyo que es más probable que al tomar un item defectuoso, provenga del turno nocturno (probabilidad aproximadamente del 60%).


###EJERCICIO 2
##Inclusión de la variable híbrido en el modelo de muestreo espacial
palmas = expand.grid(x=seq(0, 112, 7),
                     y=seq(0, 144, 9))
nrow(palmas)
## [1] 289
set.seed(123)

###Información auxiliar (producción previa)
p_racimo_u = rnorm(289, 17, 1.8)
p_racimo_p = rnorm(289, 17, 1.8)
CaMg_h17 = runif(289, 1.8, 2.0)  ##relación calcio magnesio hoja 17
CaMg_s = runif(289, 1.2, 1.4) ##relación calcio magnesio suelo
hibrid = rep(c("h1", "h2"), c(144, 145))


df= data.frame(palmas, 
               p_racimo_u, p_racimo_p,
               CaMg_h17, CaMg_s, hibrid)
head(df)
##    x y p_racimo_u p_racimo_p CaMg_h17   CaMg_s hibrid
## 1  0 0   15.99114   17.30252 1.813162 1.350780     h1
## 2  7 0   16.58568   19.10309 1.970874 1.212619     h1
## 3 14 0   19.80567   18.89753 1.902082 1.288496     h1
## 4 21 0   17.12692   19.06147 1.842000 1.300134     h1
## 5 28 0   17.23272   15.96056 1.875169 1.254049     h1
## 6 35 0   20.08712   20.60447 1.841473 1.319760     h1
###Tabla con híbridos muestreados:
  
library(clhs)
ifelse("hibrid" >= 144, "h1")
## [1] "h1"
res <- clhs(df, size=15, progress=FALSE, simple=TRUE)
## Warning: NAs introducidos por coerción
df[res,"muestreo"] = "h1"
df$muestreo[is.na(df$hibrid)] = "h2"


  
  
library(ggplot2)
 ggplot(df)+
   aes(x,y, fill=hibrid)+
   geom_tile(color="white")

##Demostrar que el modelo se acerca bastante bien a los "datos reales"  
mean(df$p_racimo_u)
## [1] 17.02191
mean(df$p_racimo_p)
## [1] 17.09966
mean(df$p_racimo_u[df$muestreo == "si"])
## [1] NA
mean(df$p_racimo_p[df$muestreo == "si"])
## [1] NA