1. Datos simulados

# Niveles de la variable x
xi <- log(seq(2500 ,4300, 200)) 
# Valores de la variable respuesta
yi <- c(10, 17, 30, 21, 18, 43, 54, 33, 60, 51) 
# Totales para cada valor de la respuesta
pesos <- c(50, 70, 100, 60, 40, 85, 90, 50, 80, 65) 
fit <- glm(yi/pesos ~ xi, weights = pesos, family = binomial(logit))
fit2 <- glm(yi/pesos ~ xi, weights = pesos, family = binomial(probit))
fit3 <- glm(yi/pesos ~ xi, weights = pesos, family = binomial(cloglog))
\(g(\mu_i)\) \(\beta_0\) \(\beta_1\) \(AIC\) \(BIC\) \(R^2_{adj}\)
\(\log (\pi /[1-\pi])\) -42.1105 5.1772 50.0784 50.6836 0.9864
\(\Phi^{-1}\left(\mu_i\right)\) -25.7556 3.1668 50.233 50.8388 0.9848
\(\log \left(-\log \left(1-\mu_i\right)\right)\) -31.0212 3.7626 49.1216 49.7268 0.9959

2. Cálculo de los límites de control

UCL <- function(be, xi, pesos, k, w = 1000){
  
  # Establecer semilla 
  set.seed(1)
  # Matriz de diseño
  X <- cbind(rep(1, length(xi)), xi)
  # Crear los objetos para cada una de las cartas
  T2H <- T2R <- T2I <- T2E <- T2D <- rep(NA, w) 
  # Vector que guardará los límites superiores
  ucls <- rep(NA, 5)
  # Matriz para los betas
  Betas <- matrix(NA, k, dim(X)[2]) 
  # Matriz para la varianza intra perfil
  S_I <- matrix(0, dim(X)[2], dim(X)[2]) 
  # Parámetro pi de la binomial
  pi_i <- expit(be[1] + be[2]*xi)
  
  # Recorrer las iteraciones
  for(j in 1:w){
    Betas[] <- NA
    S_I[] <- 0
    for(i in 1:k){
      y <- rbinom(length(xi), pesos, pi_i)
      fit <- glm(y/pesos~xi ,weights = pesos, family=binomial(logit))
      Betas[i,] <- fit$coefficients
      S_I <- S_I + vcov(fit)/k
    }
    # Guardar la media de los betas
    Bm <- colMeans(Betas)
    # Diferencia entre el parámetro estimado y su media calculada
    dif <- t(t(Betas) - Bm)
    
    # - - - - Construcción de cada carta - - - -
    
    # T2 con media muestal y matriz de covarianzas
    S_H    <- cov(Betas)
    T2Ht   <- rowSums(dif %*% chol2inv(chol(S_H)) * dif)
    T2H[j] <- max(T2Ht)
    
    # T2 con media muestral y rangos móviles
    dif2   <- Betas[2:k,] - Betas[1:(k-1),]
    S_R    <- crossprod(dif2)/(2*k-2)
    T2Rt   <- rowSums(dif %*% chol2inv(chol(S_R)) * dif)
    T2R[j] <- max(T2Rt)
    
    # T2 Intraperfil
    T2It   <- rowSums(dif %*% chol2inv(chol(S_I)) * dif)
    T2I[j] <- max(T2It)
    
    # T2 basada en MVE
    MVE    <- CovMve(Betas, nsamp = 500)
    dif3   <- t(t(Betas) - MVE$center)
    S_E    <- MVE$cov
    T2Et   <-rowSums(dif3 %*% chol2inv(chol(S_E)) * dif3)
    T2E[j] <-max(T2Et)
     
    # T2 basada en MCD
    MCD    <- cov.mcd(Betas)
    dif4   <- t(t(Betas) - MCD$center)
    S_D    <- MCD$cov
    T2Dt   <- rowSums(dif4 %*% chol2inv(chol(S_D)) * dif4)
    T2D[j] <- max(T2Dt)
  
  }
  # Guardar los limites de control superiores
  ucls[1] <- quantile(T2H, 0.95)
  ucls[2] <- quantile(T2R, 0.95)
  ucls[3] <- quantile(T2I, 0.95)
  ucls[4] <- quantile(T2E, 0.95)
  ucls[5] <- quantile(T2D, 0.95)
  return(ucls)
}

2.1. Límites encontrados tras 10000 simulaciones

\(T^2_H\) \(T^2_R\) \(T^2_I\) \(T^2_E\) \(T^2_D\)
10.73796 12.47560 12.71252 25.19752 41.06443

3. Escenarios fuera de control

3.1. Outliers

Se generaron \(k=30\) muestras, de las cuales 25 provenían del proceso bajo control y las 5 restantes de un proceso fuera de control donde el vector de parámetros del proceso fue \(\tilde{\boldsymbol{\beta}}=\boldsymbol{\beta}_0+\Delta\) siendo \(\Delta=(0.12\sigma_1,0)^{\mathrm{T}}\).

outliers <- function(be, xi, pesos, vector, k, delta_0){
  
  # Establecer semilla 
  set.seed(1)
  
  # Matriz de diseño
  X <- cbind(rep(1, length(xi)), xi)
  # Vector para los limites de control superior
  ucls <- rep(NA, 5)
  # Matriz de betas
  Betas <- matrix(NA, k, dim(X)[2]) 
  # Matriz de varianza intraperfil
  S_I <- matrix(0, dim(X)[2], dim(X)[2])
  # Parámetro pi de la binomial
  pi_i <- matrix(NA, k, length(xi))
  
  for(z in 1:k){
    pi_i[z,] <- ifelse(z%in%vector,
                       expit((be[1] + delta_0) + be[2]*xi), expit(be[1]+be[2]*xi))
  }
  # Simulación
  for(i in 1:k){
      y <- rbinom(length(xi), pesos, pi_i[i,])
      fit <- glm(y/pesos~xi, weights = pesos, family = binomial(logit))
      Betas[i,] <- fit$coefficients
      S_I <- S_I + vcov(fit)/k
  }
  # Guardar la media de los betas
  Bm <- colMeans(Betas)
  # Diferencia entre el parámetro estimado y su media calculada
  dif <- t(t(Betas) - Bm)
  
  # - - - - Construcción de cada carta - - - -
    
  # T2 con media muestal y matriz de covarianzas
  S_H <- cov(Betas)
  T2H <- rowSums(dif %*% chol2inv(chol(S_H)) * dif)
  
  # T2 con media muestral y rangos móviles
  dif2 <- Betas[2:k,] - Betas[1:(k-1),]
  S_R  <- crossprod(dif2)/(2*k-2)
  T2R  <- rowSums(dif %*% chol2inv(chol(S_R)) * dif)
  
  # T2 Intraperfil
  T2I <- rowSums(dif %*% chol2inv(chol(S_I)) * dif)
  
  # T2 basada en MVE
  MVE  <- CovMve(Betas, nsamp = 500)
  dif3 <- t(t(Betas)-MVE$center)
  S_E  <- MVE$cov
  T2E  <- rowSums(dif3 %*% chol2inv(chol(S_E)) * dif3)
  
  # T2 basada en MCD
  MCD  <- cov.mcd(Betas)
  dif4 <- t(t(Betas)-MCD$center)
  S_D  <- MCD$cov
  T2D  <- rowSums(dif4 %*% chol2inv(chol(S_D)) * dif4)
  
  resultados <- list("T2H" = T2H, "T2R" = T2R, "T2I" = T2I, "T2E" = T2E, "T2D" = T2D)
  
}

3.2. Step shifts

Se generaron \(k=30\) muestras, de las cuales la segunda mitad provenía de un proceso fuera de control donde el vector de parámetros del proceso fue \(\tilde{\boldsymbol{\beta}}=\boldsymbol{\beta}_0+\Delta\) siendo \(\Delta=(0.12\sigma_1,0)^{\mathrm{T}}\).

step_shift <- function(be, xi, pesos, k, delta_0){
  
  # Establecer semilla 
  set.seed(1)
  
  # Matriz de diseño
  X <- cbind(rep(1,length(xi)), xi)
  # Vector que guardará los límites superiores
  ucls <- rep(NA, 5)
  # Betas
  Betas <- matrix(NA, k, dim(X)[2]) 
  # Varianza intra perfil
  S_I <- matrix(0, dim(X)[2], dim(X)[2]) 
  # Parametro pi de la binomial
  pi_i <- matrix(NA, k, length(xi))
  
  for(z in 1:k){
    pi_i[z,] <- ifelse(z>k/2,
                       expit((be[1]+delta_0) + be[2]*xi), expit(be[1] + be[2]*xi))
  }
  # Simulación
  for(i in 1:k){
      y <- rbinom(length(xi), pesos, pi_i[i,])
      fit <- glm(y/pesos~xi, weights = pesos, family = binomial(logit))
      Betas[i,] <- fit$coefficients
      S_I <- S_I + vcov(fit)/k
  }
  
  # Guardar la media de los betas
  Bm <- colMeans(Betas)
  dif <- t(t(Betas) - Bm)
  
  # - - - - Construcción de cada carta - - - -
    
  # T2 con media muestal y matriz de covarianzas
  S_H <- cov(Betas)
  T2H <- rowSums(dif %*% chol2inv(chol(S_H)) * dif)
  
  # T2 con media muestral y rangos móviles
  dif2 <- Betas[2:k,] - Betas[1:(k-1),]
  S_R  <- crossprod(dif2)/(2*k-2)
  T2R  <- rowSums(dif %*% chol2inv(chol(S_R)) * dif)
  
  # T2 Intraperfil
  T2I <- rowSums(dif %*% chol2inv(chol(S_I)) * dif)
  
  # T2 basada en MVE
  MVE  <- CovMve(Betas, nsamp = 500)
  dif3 <- t(t(Betas) - MVE$center)
  S_E  <- MVE$cov
  T2E  <- rowSums(dif3 %*% chol2inv(chol(S_E)) * dif3)
  
  # T2 basada en MCD
  MCD  <- cov.mcd(Betas)
  dif4 <- t(t(Betas) - MCD$center)
  S_D  <- MCD$cov
  T2D  <- rowSums(dif4 %*% chol2inv(chol(S_D)) * dif4)
  
  resultados <- list("T2H" = T2H, "T2R" = T2R, "T2I" = T2I, "T2E" = T2E, "T2D" = T2D)

}

3.3. Drift shifts

Se generaron \(k=30\) muestras, de las cuales la segunda mitad provenía de un proceso fuera de control donde el vector de parámetros del proceso fue \(\tilde{\boldsymbol{\beta}}=\boldsymbol{\beta}_0+\frac{t-1}{k-1}\Delta\) para \(t=k/2+1,\ldots,k\), siendo \(\Delta=(0.12\sigma_1,0)^{\mathrm{T}}\).

drift_shift <- function(be, xi, pesos, k, delta_0){
  
  # Establecer semilla 
  set.seed(1)
  # Matriz de diseño
  X <- cbind(rep(1, length(xi)), xi)
  ucls <- rep(NA, 5)
  # Betas
  Betas <- matrix(NA, k, dim(X)[2]) 
  # Varianza intra perfil
  S_I <- matrix(0, dim(X)[2], dim(X)[2]) 
  # Parametro pi de la binomial
  pi_i <- matrix(NA, k, length(xi))
  
  for(z in 1:k){
    pi_i[z,] <- ifelse(z>k/2,
                       expit((be[1]+delta_0*(z-1)/(k-1)) + be[2]*xi), expit(be[1] + be[2]*xi))
  }
  for(i in 1:k){
      y <- rbinom(length(xi), pesos, pi_i[i,])
      fit <- glm(y/pesos~xi, weights = pesos, family = binomial(logit))
      Betas[i,] <- fit$coefficients
      S_I <- S_I + vcov(fit)/k
    }
  Bm <- colMeans(Betas)
  dif <- t(t(Betas) - Bm)
  
  # - - - - Construcción de cada carta - - - -
    
  # T2 con media muestal y matriz de covarianzas
  S_H <- cov(Betas)
  T2H <- rowSums(dif %*% chol2inv(chol(S_H)) * dif)
  
  # T2 con media muestral y rangos móviles
  dif2 <- Betas[2:k,] - Betas[1:(k-1),]
  S_R  <- crossprod(dif2)/(2*k-2)
  T2R  <- rowSums(dif %*% chol2inv(chol(S_R)) * dif)
  
  # T2 Intraperfil
  T2I <- rowSums(dif %*% chol2inv(chol(S_I)) * dif)
 
  # T2 basada en MVE
  MVE  <- CovMve(Betas, nsamp = 500)
  dif3 <- t(t(Betas)-MVE$center)
  S_E  <- MVE$cov
  T2E  <- rowSums(dif3 %*% chol2inv(chol(S_E)) * dif3)
  
  # T2 basada en MCD
  MCD  <- cov.mcd(Betas)
  dif4 <- t(t(Betas)-MCD$center)
  S_D  <- MCD$cov
  T2D  <- rowSums(dif4 %*% chol2inv(chol(S_D)) * dif4)
  
  resultados <- list("T2H" = T2H, "T2R" = T2R, "T2I "= T2I, "T2E" = T2E, "T2D" = T2D)

}