[R-es] Model.matrix pregunta correcta

[Javier Marcuzzi]

Estimados

En el correo anterior copie y pegue mal. Este es el correcto.

Copio y pego un código reproducible donde está la pregunta, yo lo realicé con en un archivo cmd pero si no me equivoco en rstudio o cualquier otro debería ser leído sin problemas al cortar y pegar.

---
title: "Pregunta"
author: "Javier Marcuzzi"
format: html
editor: visual
---

## Consulta y descripción del problema

De un libro estoy siguiendo unos problemas, este libro se está convirtiendo en referencia en un área muy específica, tal es así que el año pasado realizaron la cuarta edición del mismo.

Como se va convirtiendo en referencia y ser un libro donde no utilizan programas o lenguajes informáticos, distintas personas realizan en diversos lenguajes replicas de los ejemplos del libro. Tal es así, que la misma editorial publica en R la respuestas aportada por terceros.

En R tengo tres versiónes, la mía que da una falla, la de un tercero citada por la editorial y la de un Phd de una universidad de Estados Unidos.

La diferencia entre el Phd y el recomendado por la editorial, es que el primero utiliza una función para diseñar la matriz, función casi idéntica a otra universidad de Canadá, en cambio el recomendado utiliza model.matrix. Ambas formas son idénticas pero luego difieren, y esta diferencia la colocaré en ejemplos en esta consulta a la lista.

El libro tiene todo escrito en matemática, álgebra, me refiero a que no tiene ni una línea de código, por lo que cada persona es libre de utilizar papel y lápiz o la computadora.

El ejemplo en R con la función escrita por el phd de Estados Unidos, coincide con el libro, es fiel a este. Sin embargo continuando con los pasos para la resolución del problema, hay error, el mensaje es el siguiente: Error in solve.default(LHS, RHS) : Lapack routine dgesv: system is exactly singular: U\[21,21\] = 0

En cambio el otro ejemplo recomendado por la editorial, no da este error y llega a la resolución, pero no coincide con el libro, no es una copia fiel en el "armado" de la matriz.

La mejor forma de explicar la pregunta es compartiendo un código reproducible.

## Datos

Los datos son los mismos, los copio antes de continuar el trabajo, uno se llama data y el otro data0. La copia es porque al colocar "factor" no puedo utilizarlos en ambos casos. Se comprenderá en el ejemplo.

```{r}
cow = rep(seq(4, 8), each = 2)
fy = c(201, 280, 150, 200, 160, 190, 180, 250, 285, 300)
hys = c(1, 3, 1, 4, 2, 3, 1, 3, 2, 4)
pairity = rep(c(1,2), 5)

data = data.frame(cow, pairity, hys, fy)
data0 = data.frame(cow, pairity, hys, fy)
```

## Preparación de datos

Para utilizar con model.matrix, utilizo factor. En este caso se llaman data.

```{r}
data$pairity = as.factor(pairity)
data$hys = as.factor(hys)
data$hys = relevel(factor(hys), ref = "1")
data$cow = factor(data$cow)
```

Para utilizar la función de diseño, uso los datos originales (sin factor). Como no tienen "trabajo extra" para preprar los datos se llaman datos0.

## La función de diseño de matriz es la siguiente:

```{r}
# función para crear el diseño de matriz 
desgn <- function(v) {
  if (is.numeric(v)) {
    va = v
    mrow = length(va)
    mcol = max(va)
  }
  if (is.character(v)) {
    vf = factor(v)
    # Guarda el indice para cada nivel
    va = as.numeric(vf)
    mrow = length(va)
    mcol = length(levels(vf))
  }
  
  # Inicio la martiz
  X = matrix(data = c(0), nrow = mrow, ncol = mcol)
  
  for (i in 1:mrow) {
    ic = va[i]
    X[i, ic] = 1
  }
  return(X)
}
```

## Armado de las matrices

### Forma 1

Utilizando la función de diseño, hace una fiel coincidencia con el libro, me refiero a que el libro tiene todo escrito en forma matemática, sin utilizar lenguajes informáticos, se puede hacer con lápiz y papel, no hay otro requerimiento. Pero luego este armado de matrices da error (singularidad).

Como se puede reproducir, el ejemplo son dos matrices que se unen con cbind.

```{r}
# Matriz X, relaciona los registros de HYS y partos
# dideño para HYS
x1 = desgn(data0$hys)
x1
# diseño para Partos
x2 = desgn(data0$pairity)
x2
# Matriz X segín libro
XLibro = cbind(x1,x2)

t(XLibro)
```

### Forma 2

Utilizando los datos con model.matrix, de esta forma también luego el mismo error (singularidad) y no coincide con el libro.

```{r}
X = model.matrix(fy ~ -1 + pairity + hys, data = data)
X
```

### Forma 3 X\[,-4\]

La diferencia con el anterior está en X\[,-4\], como este es el que funciona lo llamo XOK.

Funciona quiere decir que no da error, los resultados coinciden con el libro, pero no visualmente cuando leo las matrices. Este punto se ejemplifica en la comparación. La única diferencia con los resultados del libro, es que "elimina" los dos casos donde da 0,000.

```{r}
XOK = X[,-4]
```

## Comparación

Para analizar comparo las tres formas de armado de la matriz, las dos primeras luego dan singularidad y la tercera no,

### Libro pero que da singularidad (forma 1)

```{r}
XLibro
```

### model.matrix que también da singularidad (forma 2)

```{r}
X
```

### model.matrix que da los resultados correctos (forma 3 X\[,-4\])

```{r}
XOK
```

## Pregunta

Si se observa, Xlibro, X y XOK, todos tienen 10 filas, esto está bién, pero Xlibro tiene 6 columnas, X tiene 5 columnas y XOK tiene 4 columnas.

El de 6 columnas es identico al libro, pero da singularidad, el de 4 columnas da los resultados del libro, con una pequeña diferencia en los dos que dan 0,000, pero el "armado" con model.matrix difiere.

¿Alguien puede explicarme model.matrix? En el diseño, o en la forma de armar optimiza algo, pero no comprendo como realiza la convinación, en forma automática eliminina \[,-4\] los dos que darían 0.

Recuerdo, la formula es fy \~ -1 + pairity + hys o dos matrices unidas con cbind (libro).

¿Cómo leer lo siguiente?

```{r}
XLibro
X
XOK
```

##