Análisis de Datos 2 - examen extraordinario - 2020-jun-30

Justifica las respuestas.

1. (1 punto) Sea \(X\) una variable aleatoria con distribución de Poisson P(\(\lambda=5\)).

   Sea \(Y\) una variable aleatoria tal que \(Y=X+\theta\), donde \(\theta\) es un número entero.

   Obtén la estimación máximo-verosímil de \(\theta\) para la muestra de \(Y\) siguiente:

   \[ y_1=6 \quad y_2=7 \quad y_3=7 \quad y_4=12 \quad y_5=3 \]

   vero1 <- function (x, zita) dpois (x-zita, 5)
   vero  <- function (zita) prod (sapply (c(6,7,7,12,3), function(x) vero1(x,zita)))
   posibles <- (-10):10
   which.max (structure (sapply (posibles, vero), names=posibles))
   # la estimación MV de zita es 2
   

2. (1 punto) Se pretende servir una comida para trescientas personas minimizando los gastos. Cada persona requiere consumir 200 kcal, 30 g proteína, 1 g lípido y 1 g glúcido. La siguiente tabla presenta cantidades por cada 100 g de alimento:

   alimento	kcal	g proteína	g lípido	g glúcido	€
   lácteo	60	4	1,8	3,6	1,7
   carne	180	22	5,0	0,4	2,5
   pescado	82	17	1,3	1,2	3,1
   legumbre	100	15	3,4	3,3	1,8
   fruta	48	2	2,7	9,8	0,8

   ¿Cuánto hay que adquirir de cada alimento?

   library (lpSolve)
   sol <- lp (objective = c(1.7,2.5,3.1,1.8,0.8),
              const.mat = t (matrix (c (60,180,82,100,48,
                                        4,22,17,15,2,
                                        1.8,5.0,1.3,3.4,2.7,
                                        3.6,0.4,1.2,3.3,9.8),
                                     5)),
              const.dir = ">",
              const.rhs = 300 * c(200,30,1,1))
   structure (0.1*sol$solution, names=c("lac","car","pes","leg","fru"))
   #                        lac       car       pes       leg       fru 
   # en quilogramos:   0.000000 37.837838  0.000000  4.504505  0.000000
   

3. (3 puntos) Plantea y resuelve el problema anterior como un recocido simulado.

   ### usando optim 
   A <- t (matrix (c (60,180,82,100,48,
                      4,22,17,15,2,
                      1.8,5.0,1.3,3.4,2.7,
                      3.6,0.4,1.2,3.3,9.8),
                   5))
   b <- 300 * c(200,30,1,1)
   fun <- function (x, penal=10000)
   {
      lac <- x[1]; car <- x[2]; pes <- x[3]; leg <- x[4]; fru <- x[5]
      1.7 * lac + 2.5 * car + 3.1 * pes + 1.8 * leg + 0.8 * fru +
        sum (A %*% x < b) * penal +
        sum (x < 0) * penal
   }
   optim (rep(1000,5), fun, method="SANN", control = list(trace=1,maxit=1e6))
   # [1]   0.07240620 236.09558366   0.21216099 253.50925224   0.01873769
   # $value  [1] 1047.351
   optim (rep(1000,5), fun, method="SANN", control = list(trace=1,maxit=1e8,temp=1000))
   # [1]   1.570634 338.577090   2.497258  99.619832   4.846377
   # $value  [1] 1040.047
   
   ### usando programa propio
   f <- function (lac, car, pes, leg, fru, penal = 10000) 
          1.7 * lac + 2.5 * car + 3.1 * pes + 1.8 * leg + 0.8 * fru +
          sum (A %*% c(lac,car,pes,leg,fru) < b) * penal +
          sum (c(lac,car,pes,leg,fru) < 0) * penal
   xmin <- 0; xmax <- 1000
   T <- 10; sigma <- 0.5
   lac0 <- car0 <- pes0 <- leg0 <- fru0 <- 1000
   y0 <- f(lac0,car0,pes0,leg0,fru0)
   for (n in 1:1000000) 
   {
    T <- T * 0.9999
    x1 <- pmax (xmin, pmin (xmax, rnorm (5, c(lac0,car0,pes0,leg0,fru0), sigma)))
    y1 <- do.call (f, as.list(x1))
    if (y1 <= y0) alfa <- 1 else
      alfa <- exp ((y0-y1)/T)
    if (runif(1) < alfa) 
      {
         lac0 <- x1[1]
         car0 <- x1[2]
         pes0 <- x1[3]
         leg0 <- x1[4]
         fru0 <- x1[5]
         y0 <- y1
      }
   }
   c(lac=lac0,car=car0,pes=pes0,leg=leg0,fru=fru0,y=y0)
   ##     lac        car        pes        leg        fru          y 
   ## 0.00000  378.37615    0.00000   45.04985    0.00000 1027.03010
   

4. (2 puntos) Paraleliza la implementación del problema anterior.

   f <- function (lac, car, pes, leg, fru, penal = 10000) 
          1.7 * lac + 2.5 * car + 3.1 * pes + 1.8 * leg + 0.8 * fru +
          sum (A %*% c(lac,car,pes,leg,fru) < b) * penal +
          sum (c(lac,car,pes,leg,fru) < 0) * penal
   xmin <- 0; xmax <- 1000
   T <- 10; sigma <- 0.5
   lac0 <- car0 <- pes0 <- leg0 <- fru0 <- 1000
   y0 <- f(lac0,car0,pes0,leg0,fru0)
   unrecocido <- function (i) # i no se usará
   {
      for (n in 1:1000000) 
      {
       T <- T * 0.9999
       x1 <- pmax (xmin, pmin (xmax, rnorm (5, c(lac0,car0,pes0,leg0,fru0), sigma)))
       y1 <- do.call (f, as.list(x1))
       if (y1 <= y0) alfa <- 1 else
         alfa <- exp ((y0-y1)/T)
       if (runif(1) < alfa) 
         {
            lac0 <- x1[1]
            car0 <- x1[2]
            pes0 <- x1[3]
            leg0 <- x1[4]
            fru0 <- x1[5]
            y0 <- y1
         }
      }
      c(lac=lac0,car=car0,pes=pes0,leg=leg0,fru=fru0,y=y0)
   }
   library (parallel)
   nhilos <- detectCores()-1
   soluciones <- mclapply(1:nhilos, unrecocido, mc.cores=nhilos)
   soluciones [[which.min (sapply (soluciones, function(s) s["y"]))]]