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:

###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)

###Cambio el tamaño a 20 porque con 15 habían cuadros sin muestrear
res = clhs(df, size=20, progress = FALSE, simple = TRUE)
## Warning: NAs introducidos por coerción
resH = ifelse(res>= 144, "h1", "h2")

df[res, 'muestreo'] = resH
df$muestreo[is.na(df$muestreo)] = 'sin muestrear'

table(resH)
## resH
## h1 h2 
## 11  9
library(ggplot2)

ggplot(df)+
  aes(x,y, fill= muestreo)+
  geom_tile(color = 'white')

###Cálculo de medias:
mean(df$p_racimo_u)
## [1] 17.02191
mean(df$p_racimo_p)
## [1] 17.09966
mean(df$p_racimo_u[df$muestreo == "h1"])
## [1] 17.2233
mean(df$CaMg_h17[df$muestreo == "h1"])
## [1] 1.909184
mean(df$p_racimo_u[df$muestreo == "h2"])
## [1] 17.06161
mean(df$CaMg_h17[df$muestreo == "h2"])
## [1] 1.891654
##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 == "h1"])
## [1] 17.2233
mean(df$p_racimo_p[df$muestreo == "h2"])
## [1] 17.48285