Cómo realizar la transformación logarítmica isométrica

Question

Tengo datos sobre los comportamientos de movimiento (tiempo dedicado a dormir, al sedentarismo y a la actividad física) que suman aproximadamente 24 (en horas al día). Quiero crear una variable que capte el tiempo relativo dedicado a cada uno de estos comportamientos; me han dicho que una transformación isométrica de logaritmo lo lograría.

Parece que debería usar la función ilr en R, pero no encuentro ningún ejemplo real con código. ¿Por dónde empiezo?

Las variables que tengo son el tiempo de sueño, el promedio de tiempo sedentario, el promedio de actividad física ligera, el promedio de actividad física moderada y el promedio de actividad física vigorosa. El tiempo de sueño fue autodeclarado, mientras que las otras son promedios de días válidos de datos del acelerómetro. Así que para estas variables, los casos no suman exactamente 24.

Mi suposición: Estoy trabajando en SAS, pero parece que R será mucho más fácil de usar para esta parte. Así que primero importar los datos con sólo las variables de interés. Luego usar la función acomp(). Luego no consigo entender la sintaxis de la función ilr(). Cualquier ayuda será muy apreciada.

Answer 1

La transformación ILR (Isometric Log-Ratio) se utiliza en el análisis de los datos de composición. Cualquier observación es un conjunto de valores positivos que suman la unidad, como las proporciones de sustancias químicas en una mezcla o las proporciones del tiempo total dedicado a diversas actividades. La invariante suma-unidad implica que, aunque puede haber $k\ge 2$ componentes a cada observación, sólo hay $k-1$ valores funcionalmente independientes. (Geométricamente, las observaciones se encuentran en un $k-1$ -simplex de dimensiones en $k$ -espacio euclidiano de dimensiones $\mathbb{R}^k$ . Esta simplicidad se manifiesta en las formas triangulares de los gráficos de dispersión de los datos simulados que se muestran a continuación).

Normalmente, las distribuciones de los componentes se vuelven más "agradables" cuando se transforman en logaritmos. Esta transformación puede escalarse dividiendo todos los valores de una observación por su media geométrica antes de tomar los logaritmos. (Equivalentemente, los logaritmos de los datos de cualquier observación se centran restando su media). Esto se conoce como la transformación "Centered Log-Ratio", o CLR. Los valores resultantes siguen estando dentro de un hiperplano en $\mathbb{R}^k$ porque el escalado hace que la suma de los logaritmos sea cero. El ILR consiste en elegir cualquier base ortonormal para este hiperplano: el $k-1$ Las coordenadas de cada observación transformada se convierten en sus nuevos datos. De forma equivalente, el hiperplano se gira (o se refleja) para que coincida con el plano con fuga $k^\text{th}$ y se utiliza la primera $k-1$ coordenadas. (Como las rotaciones y las reflexiones conservan la distancia, son isometrías de ahí el nombre de este procedimiento).

Tsagris, Preston y Wood afirman que "una elección estándar de [la matriz de rotación] $H$ es la submatriz de Helmert obtenida al eliminar la primera fila de la matriz de Helmert".

La matriz de Helmert de orden $k$ se construye de manera sencilla (véase Harville p. 86, por ejemplo). Su primera fila es toda $1$ s. La siguiente fila es una de las más sencillas que puede hacerse ortogonal a la primera fila, a saber $(1, -1, 0, \ldots, 0)$ . Fila $j$ es uno de los más sencillos que es ortogonal a todas las filas precedentes: su primera $j-1$ las entradas son $1$ s, lo que garantiza que es ortogonal a las filas $2, 3, \ldots, j-1$ y su $j^\text{th}$ la entrada se establece en $1-j$ para que sea ortogonal a la primera fila (es decir, sus entradas deben sumar cero). A continuación, todas las filas se reescalan a la unidad de longitud.

Aquí, para ilustrar el patrón, está el $4\times 4$ Matriz de Helmert antes de que sus filas hayan sido reescaladas:

$$\pmatrix{1&1&1&1 \\ 1&-1&0&0 \\ 1&1&-2&0 \\ 1&1&1&-3}.$$

(Edición añadida en agosto de 2017) Un aspecto especialmente agradable de estos "contrastes" (que se leen fila a fila) es su interpretabilidad. La primera fila se elimina, dejando $k-1$ las filas restantes para representar los datos. La segunda fila es proporcional a la diferencia entre la segunda variable y la primera. La tercera fila es proporcional a la diferencia entre la tercera variable y las dos primeras. En general, la fila $j$ ( $2\le j \le k$ ) refleja la diferencia entre la variable $j$ y todas las que le preceden, las variables $1, 2, \ldots, j-1$ . Esto deja la primera variable $j=1$ como "base" para todos los contrastes. Estas interpretaciones me han resultado útiles al seguir el ILR mediante el Análisis de Componentes Principales (ACP): permite interpretar las cargas, al menos de forma aproximada, en términos de comparaciones entre las variables originales. He insertado una línea en el R aplicación de ilr que da a las variables de salida nombres adecuados para ayudar a esta interpretación. (Fin de la edición.)

Desde R proporciona una función contr.helmert para crear tales matrices (aunque sin el escalado, y con filas y columnas negadas y transpuestas), ni siquiera tienes que escribir el (simple) código para hacerlo. Usando esto, implementé el ILR (ver abajo). Para ejercitarlo y probarlo, generé $1000$ extracciones independientes de una distribución Dirichlet (con parámetros $1,2,3,4$ ) y trazamos su matriz de dispersión. Aquí, $k=4$ .

Todos los puntos se agrupan cerca de las esquinas inferiores izquierdas y llenan parches triangulares de sus áreas de trazado, como es característico de los datos de composición.

Su ILR tiene sólo tres variables, de nuevo trazadas como una matriz de dispersión:

En efecto, el aspecto es más bonito: los gráficos de dispersión han adquirido formas más características de "nube elíptica", que se adaptan mejor a los análisis de segundo orden, como la regresión lineal y el ACP.

Tsagris y otros. generalizar el CLR utilizando una transformación Box-Cox, que generaliza el logaritmo. (El logaritmo es una transformación Box-Cox con parámetro $0$ .) Es útil porque, como argumentan los autores (correctamente en mi opinión), en muchas aplicaciones los datos deberían determinar su transformación. Para estos datos de Dirichlet un parámetro de $1/2$ (que está a medio camino entre la no transformación y la transformación logarítmica) funciona de maravilla:

"Hermoso" se refiere a la simple descripción que permite esta imagen: en lugar de tener que especificar la ubicación, la forma, el tamaño y la orientación de cada nube de puntos, sólo tenemos que observar que (con una excelente aproximación) todas las nubes son circulares con radios similares. En efecto, el CLR ha simplificado una descripción inicial que requería al menos 16 números en una que sólo requiere 12 números y el ILR ha reducido eso a sólo cuatro números (tres ubicaciones univariantes y un radio), al precio de especificar el parámetro ILR de $1/2$ --un quinto número. Cuando estas simplificaciones dramáticas ocurren con datos reales, normalmente pensamos que estamos en algo: hemos hecho un descubrimiento o ha conseguido un la perspicacia.

---

Esta generalización se aplica en el ilr función de abajo. El comando para producir estas variables "Z" era simplemente

z <- ilr(x, 1/2)

Una ventaja de la transformación de Box-Cox es su aplicabilidad a las observaciones que incluyen verdaderos ceros: sigue estando definida siempre que el parámetro sea positivo.

Referencias

Michail T. Tsagris, Simon Preston y Andrew T.A. Wood, Una transformación de potencia basada en datos para los datos de composición . arXiv:1106.1451v2 [stat.ME] 16 de junio de 2011.

David A. Harville, Álgebra matricial desde la perspectiva del estadístico . Springer Science & Business Media, 27 de junio de 2008.

---

Aquí está el R código.

#
# ILR (Isometric log-ratio) transformation.
# `x` is an `n` by `k` matrix of positive observations with k >= 2.
#
ilr <- function(x, p=0) {
  y <- log(x)
  if (p != 0) y <- (exp(p * y) - 1) / p       # Box-Cox transformation
  y <- y - rowMeans(y, na.rm=TRUE)            # Recentered values
  k <- dim(y)[2]
  H <- contr.helmert(k)                       # Dimensions k by k-1
  H <- t(H) / sqrt((2:k)*(2:k-1))             # Dimensions k-1 by k
  if(!is.null(colnames(x)))                   # (Helps with interpreting output)
    colnames(z) <- paste0(colnames(x)[-1], ".ILR")
  return(y %*% t(H))                          # Rotated/reflected values
}
#
# Specify a Dirichlet(alpha) distribution for testing.
#
alpha <- c(1,2,3,4)
#
# Simulate and plot compositional data.
#
n <- 1000
k <- length(alpha)
x <- matrix(rgamma(n*k, alpha), nrow=n, byrow=TRUE)
x <- x / rowSums(x)
colnames(x) <- paste0("X.", 1:k)
pairs(x, pch=19, col="#00000040", cex=0.6)
#
# Obtain the ILR.
#
y <- ilr(x)
colnames(y) <- paste0("Y.", 1:(k-1))
#
# Plot the ILR.
#
pairs(y, pch=19, col="#00000040", cex=0.6)

Answer 2

Para su caso de uso, probablemente esté bien reducir todo a uno. El hecho de que los números no sumen exactamente 24 añadirá un poco de ruido adicional a los datos, pero no debería estropear mucho las cosas.

Como bien dijo @whuber, al tratarse de proporciones, hay que tener en cuenta las dependencias entre las variables (ya que suman una). La transformada ilr se ocupa adecuadamente de esto, ya que transforma las variables en $\mathbb{R}^{D-1}$ para $D$ proporciones.

Dejando de lado todos los detalles técnicos, es importante saber interpretar correctamente los datos transformados por la ilr. A fin de cuentas, la transformada ilr sólo se refiere a los cocientes logarítmicos de los grupos. Pero la define con respecto a alguna jerarquía predefinida. Si se define una jerarquía como la siguiente

cada variable transformada puede calcularse como

$b_i = \sqrt{\frac{rs}{r + s}}\ln \frac{g(R_i)}{g(S_i)}$

donde $i$ representa un nodo interno de la jerarquía, $R_i$ define una partición de variables correspondiente a $i$ , $S_i$ define la otra partición de variables correspondiente a $i$ y $g(...)$ se refiere a la media geométrica. Estas variables transformadas también se conocen como saldos.

La siguiente pregunta es: ¿cómo se define la jerarquía de las variables? Esto depende realmente de ti, pero si tienes tres variables, no hay demasiadas combinaciones con las que jugar. Por ejemplo, puedes definir la jerarquía como

                        /-A
            /(A|B)-----|
-(AB|C)----|            \-B
           |
            \-C

donde A representa el tiempo de sueño, B representa el tiempo de sedentarismo, C representa el tiempo dedicado a la actividad física (A|B) representa la relación logarítmica normalizada entre $A$ $B$ (es decir $\frac{1}{\sqrt{2}}\ln \frac{A}{B}$ ), y $(AB|C)$ se refiere a la relación logarítmica normalizada entre $A$ , $B$ y $C$ (es decir $\frac{\sqrt{2}}{\sqrt{3}} \ln \frac{AB}{C}$ ). Si hay muchas variables, compruebo algunos de los trabajos realizados con saldos principales

Pero volviendo a tu pregunta original, ¿cómo puedes utilizar esta información para realizar realmente la transformación ilr?

Si estás usando R, yo revisaría el paquete de composiciones

Para usar ese paquete, tendrás que entender cómo crear una partición binaria secuencial (SBP), que es como se define la jerarquía. Para la jerarquía definida anteriormente, puedes representar la SBP con la siguiente matriz.

        A  B  C
(A|B)   1 -1  0
(AB|C)  1  1 -1

donde los valores positivos representan las variables del numerador, los valores negativos representan las variables del denominador y los ceros representan la ausencia de esa variable en la balanza. Se puede construir la base ortonormal utilizando balanceBase del SBP que ha definido.
Una vez que tengas esto deberías poder pasar tu tabla de proporciones junto con la base que calculaste anteriormente.

Yo revisaría esto referencia para la definición original de los saldos

Answer 3

Las entradas anteriores responden a la pregunta sobre cómo construir un base de ILR y obtener sus saldos de ILR. Además, la elección de que base puede facilitar la interpretación de sus resultados.

Es posible que le interese una partición la siguiente partición:

(1) (dormir,sedentarismo|actividad_física) (2) (dormir|sedentario).

Como tiene tres partes en su composición, obtendrá dos balances ILR para analizar. Configurando la partición como se ha indicado anteriormente, podrá obtener los balances correspondientes a "activo o no" (1) y "qué forma de inactividad" (2).

Si se analiza cada saldo de ILR por separado, por ejemplo realizando una regresión con respecto a la hora del día o del año para ver si hay cambios, se pueden interpretar los resultados en términos de cambios en "activo o no" y cambios en "qué forma de inactividad".

Si, por el contrario, va a realizar técnicas como el PCA que obtienen una nueva base en el espacio ILR, sus resultados no dependerán de la elección de la partición. Esto se debe a que sus datos existen en el espacio CLR, el plano D-1 ortogonal al vector uno, y los equilibrios ILR son diferentes opciones de ejes de norma unitaria para describir la posición de los datos en el plano CLR.