arrays - una - size matlab para que sirve

Metodología de estructuración de datos de alta dimensión en R vs. MATLAB

Question

Pregunta

¿Cuál es la forma correcta de estructurar datos multivariantes con etiquetas categóricas acumuladas en ensayos repetidos para el análisis exploratorio en R? No quiero volver a MATLAB.

Explicación

Me gustan mucho más las funciones de análisis y la sintaxis de R (y las tramas asombrosas) que las de MATLAB, y he estado trabajando duro para refactorizar mis cosas. Sin embargo, sigo colgando sobre la forma en que se organizan los datos en mi trabajo.

MATLAB

Es típico para mí trabajar con series de tiempo multivariantes repetidas en muchos ensayos, que se almacenan en una gran matriz multidimensional rango-3 tensor de SERIESxSAMPLESxTRIALS. Esto se presta a algunas cosas buenas de álgebra lineal de vez en cuando, pero es torpe cuando se trata de otra variable, a saber, CLASE. Normalmente, las etiquetas de clase se almacenan en otro vector de dimensión 1x TRIALS .

Cuando se trata de análisis, trazo lo menos posible, porque cuesta tanto trabajo reunir una trama realmente buena que te enseña mucho sobre los datos en MATLAB. ( No soy el único que se siente de esta manera ).

R

En R, me he mantenido lo más cerca posible de la estructura de MATLAB, pero las cosas se vuelven molestamente complejas cuando trato de mantener el etiquetado de la clase por separado; Tendría que seguir pasando las etiquetas en funciones a pesar de que solo estoy usando sus atributos. Entonces lo que hice fue separar la matriz en una lista de matrices por CLASS. Esto agrega complejidad a todas mis funciones de apply() , pero parece valer la pena en términos de mantener las cosas consistentes (y errores).

Por otro lado, R simplemente no parece ser amigable con tensores / matrices multidimensionales. Solo para trabajar con ellos, necesitas agarrar la abind los abind . La documentación sobre análisis multivariante, como este ejemplo, parece operar bajo la suposición de que tienes una enorme tabla 2-D de puntos de datos, como un pergamino medieval largo, un marco de datos, y no menciona cómo llegar "allí" desde donde estoy. .

Una vez que logro trazar y clasificar los datos procesados, no es un problema tan grande, ya que para entonces he trabajado en estructuras de marcos de datos con formas como TRIALSxFEATURES (la melt ha ayudado mucho con esto). Por otro lado, si quiero generar rápidamente una matriz de diagramas de dispersión o un histograma de latticista para la fase exploratoria (es decir, momentos estadísticos, separación, varianza entre clases, histogramas, etc.), debo detenerme y descubrir cómo Voy a apply() estas enormes matrices multidimensionales en algo que entienden esas bibliotecas.

Si continúo latiendo en la jungla con soluciones ad-hoc para esto, o nunca voy a mejorar o terminaré con mis propias y extrañas maneras de hacerlo que no tienen sentido para nadie .

Entonces, ¿cuál es la forma correcta de estructurar datos multivariantes con etiquetas categóricas acumuladas en ensayos repetidos para el análisis exploratorio en R? Por favor, no quiero regresar a MATLAB.

Bonificación: tiendo a repetir estos análisis sobre estructuras de datos idénticas para múltiples sujetos. ¿Hay una mejor manera general de envolver los fragmentos de código en los bucles?

Answer 1

Tal vez dplyr :: tbl_cube?

Partiendo de la excelente respuesta de @ BrodieG, creo que puede resultarle útil observar la nueva funcionalidad disponible en dplyr::tbl_cube . Esto es esencialmente un objeto multidimensional que puedes crear fácilmente a partir de una lista de matrices (como lo estás haciendo en la actualidad), que tiene algunas funciones realmente buenas para subconjuntos, filtrado y resumen que (importante, creo) se usan de manera consistente en todo el " cubo "vista" y "tabular" vista de los datos.

require(dplyr)

Un par de advertencias:

Es un lanzamiento temprano: todos los problemas que van con eso
Se recomienda para esta versión descargar plyr cuando se carga dplyr

Cargando matrices en cubos

Aquí hay un ejemplo usando arr como se define en la otra respuesta:

# using arr from previous example
# we can convert it simply into a tbl_cube
arr.cube<-as.tbl_cube(arr)

arr.cube  
#Source: local array [24 x 3]  
#D: ser [chr, 3]  
#D: smp [chr, 2]  
#D: tr [chr, 4]  
#M: arr [dbl[3,2,4]]

Tenga en cuenta que D significa Dimensiones y M Medidas, y puede tener tantas como quiera de cada una.

Fácil conversión de multidimensional a plana

Puede convertir fácilmente la tabla de datos devolviéndola como un data.frame (que simplemente puede convertir a un data.table si necesita la funcionalidad y los beneficios de rendimiento más adelante)

head(as.data.frame(arr.cube))
#    ser   smp   tr       arr
#1 ser 1 smp 1 tr 1 0.6656456
#2 ser 2 smp 1 tr 1 0.6181301
#3 ser 3 smp 1 tr 1 0.7335676
#4 ser 1 smp 2 tr 1 0.9444435
#5 ser 2 smp 2 tr 1 0.8977054
#6 ser 3 smp 2 tr 1 0.9361929

Subconjunto

Obviamente, podría aplanar todos los datos para cada operación, pero eso tiene muchas implicaciones para el rendimiento y la utilidad. Creo que el beneficio real de este paquete es que puede "preexaminar" el cubo para los datos que necesita antes de convertirlo en un formato tabular que sea ggplot-friendly, por ejemplo, un filtrado simple para devolver solo la serie 1:

arr.cube.filtered<-filter(arr.cube,ser=="ser 1")
as.data.frame(arr.cube.filtered)
#    ser   smp   tr       arr
#1 ser 1 smp 1 tr 1 0.6656456
#2 ser 1 smp 2 tr 1 0.9444435
#3 ser 1 smp 1 tr 2 0.4331116
#4 ser 1 smp 2 tr 2 0.3916376
#5 ser 1 smp 1 tr 3 0.4669228
#6 ser 1 smp 2 tr 3 0.8942300
#7 ser 1 smp 1 tr 4 0.2054326
#8 ser 1 smp 2 tr 4 0.1006973

tbl_cube actualmente trabaja con las funciones dplyr summarise() , select() , group_by() y filter() . De manera útil, puede encadenarlos junto con el operador %.% .

Para el resto de los ejemplos, voy a usar el objeto incorporado nasa tbl_cube, que tiene un montón de datos meteorológicos (y muestra múltiples dimensiones y medidas):

Agrupación y medidas sumarias

nasa
#Source: local array [41,472 x 4]
#D: lat [dbl, 24]
#D: long [dbl, 24]
#D: month [int, 12]
#D: year [int, 6]
#M: cloudhigh [dbl[24,24,12,6]]
#M: cloudlow [dbl[24,24,12,6]]
#M: cloudmid [dbl[24,24,12,6]]
#M: ozone [dbl[24,24,12,6]]
#M: pressure [dbl[24,24,12,6]]
#M: surftemp [dbl[24,24,12,6]]
#M: temperature [dbl[24,24,12,6]]

Así que aquí hay un ejemplo que muestra lo fácil que es extraer un subconjunto de datos modificados del cubo, y luego aplanarlo para que sea apropiado para el trazado:

plot_data<-as.data.frame(          # as.data.frame so we can see the data
filter(nasa,long<(-70)) %.%        # filter long < (-70) (arbitrary!)
group_by(lat,long) %.%             # group by lat/long combo
summarise(p.max=max(pressure),     # create summary measures for each group
          o.avg=mean(ozone),
          c.all=(cloudhigh+cloudlow+cloudmid)/3)
)

head(plot_data)

#       lat   long p.max    o.avg    c.all
#1 36.20000 -113.8   975 310.7778 22.66667
#2 33.70435 -113.8   975 307.0833 21.33333
#3 31.20870 -113.8   990 300.3056 19.50000
#4 28.71304 -113.8  1000 290.3056 16.00000
#5 26.21739 -113.8  1000 282.4167 14.66667
#6 23.72174 -113.8  1000 275.6111 15.83333

Notación consistente para estructuras de datos nd y 2-d

Lamentablemente, la función tbl_cube mutate() aún no está implementada para tbl_cube pero parece que será cuestión de (no mucho) tiempo. Puede usarlo (y todas las otras funciones que funcionan en el cubo) en el resultado tabular, sin embargo, con exactamente la misma notación. Por ejemplo:

plot_data.mod<-filter(plot_data,lat>25) %.%    # filter out lat <=25
mutate(arb.meas=o.avg/p.max)                   # make a new column

head(plot_data.mod)

#       lat      long p.max    o.avg    c.all  arb.meas
#1 36.20000 -113.8000   975 310.7778 22.66667 0.3187464
#2 33.70435 -113.8000   975 307.0833 21.33333 0.3149573
#3 31.20870 -113.8000   990 300.3056 19.50000 0.3033389
#4 28.71304 -113.8000  1000 290.3056 16.00000 0.2903056
#5 26.21739 -113.8000  1000 282.4167 14.66667 0.2824167
#6 36.20000 -111.2957   930 313.9722 20.66667 0.3376045

Trazado: como un ejemplo de la funcionalidad R que le gusta a los datos planos

Luego puede trazar con ggplot() usando los beneficios de los datos aplanados:

# plot as you like:
ggplot(plot_data.mod) +
  geom_point(aes(lat,long,size=c.all,color=c.all,shape=cut(p.max,6))) +
  facet_grid( lat ~ long ) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1))

Usando data.table en los datos planos resultantes

No voy a ampliar el uso de data.table aquí, ya que está bien en la respuesta anterior. Obviamente hay muchas buenas razones para usar data.table : para cualquier situación aquí puede devolver una mediante una simple conversión de data.frame:

data.table(as.data.frame(your_cube_name))

Trabajando dinámicamente con tu cubo

Otra cosa que creo que es genial es la capacidad de agregar medidas (cortes / escenarios / cambios, como quiera llamarlas) a su cubo. Creo que esto encajará bien con el método de análisis descrito en la pregunta. Aquí hay un ejemplo simple con arr.cube , que agrega una medida adicional que es en sí misma una función (ciertamente simple) de la medida anterior. Accede / actualiza medidas a través de la sintaxis yourcube $mets[$...]

head(as.data.frame(arr.cube))

#    ser   smp   tr       arr
#1 ser 1 smp 1 tr 1 0.6656456
#2 ser 2 smp 1 tr 1 0.6181301
#3 ser 3 smp 1 tr 1 0.7335676
#4 ser 1 smp 2 tr 1 0.9444435
#5 ser 2 smp 2 tr 1 0.8977054
#6 ser 3 smp 2 tr 1 0.9361929

arr.cube$mets$arr.bump<-arr.cube$mets$arr*1.1  #arb modification!

head(as.data.frame(arr.cube))

#    ser   smp   tr       arr  arr.bump
#1 ser 1 smp 1 tr 1 0.6656456 0.7322102
#2 ser 2 smp 1 tr 1 0.6181301 0.6799431
#3 ser 3 smp 1 tr 1 0.7335676 0.8069244
#4 ser 1 smp 2 tr 1 0.9444435 1.0388878
#5 ser 2 smp 2 tr 1 0.8977054 0.9874759
#6 ser 3 smp 2 tr 1 0.9361929 1.0298122

Dimensiones - o no ...

He jugado un poco al tratar de agregar dinámicamente dimensiones completamente nuevas (ampliando efectivamente un cubo existente con dimensiones adicionales y clonando o modificando los datos originales usando yourcube $dims[$...] dims $dims[$...] ) pero he encontrado que el comportamiento es un poco inconsistente. Probablemente sea mejor evitar esto de todos modos, y estructure su cubo primero antes de manipularlo. Te mantendré informado si llego a algún lado.

Persistencia

Obviamente, uno de los principales problemas para tener acceso de intérprete a una base de datos multidimensional es la posibilidad de incomodarlo accidentalmente con una pulsación de tecla inoportuna. Así que supongo que persistir temprano y con frecuencia:

tempfilename<-gsub("[ :-]","",paste0("DBX",(Sys.time()),".cub"))
# save:
save(arr.cube,file=tempfilename)
# load:
load(file=tempfilename)

¡Espero que ayude!

Answer 2

Como se ha señalado, muchas de las herramientas analíticas y de visualización más potentes se basan en datos en formato largo. Ciertamente, para las transformaciones que se benefician del álgebra matricial, debe mantener las cosas en matrices, pero tan pronto como desee ejecutar análisis paralelos en subconjuntos de sus datos, o trazar cosas por factores en sus datos, realmente desea melt .

Aquí hay un ejemplo para comenzar con data.table y ggplot .

Matriz -> Tabla de datos

Primero, hagamos algunos datos en su formato:

series <- 3
samples <- 2
trials <- 4

trial.labs <- paste("tr", seq(len=trials))
trial.class <- sample(c("A", "B"), trials, rep=T)

arr <- array(
  runif(series * samples * trials), 
  dim=c(series, samples, trials),
  dimnames=list(
    ser=paste("ser", seq(len=series)), 
    smp=paste("smp", seq(len=samples)), 
    tr=trial.labs
  )
)
# , , tr = Trial 1
#        smp
# ser         smp 1     smp 2
#   ser 1 0.9648542 0.4134501
#   ser 2 0.7285704 0.1393077
#   ser 3 0.3142587 0.1012979
#
# ... omitted 2 trials ...
# 
# , , tr = Trial 4
#        smp
# ser         smp 1     smp 2
#   ser 1 0.5867905 0.5160964
#   ser 2 0.2432201 0.7702306
#   ser 3 0.2671743 0.8568685

Ahora tenemos una matriz tridimensional. data.table y data.table en un data.table (note que melt en data.frames , que son básicamente data.table s sans bells & whistles, entonces primero tenemos que fundir, luego convertir a data.table ):

library(reshape2)
library(data.table)

dt.raw <- data.table(melt(arr), key="tr")  # we''ll get to what the `key` arg is doing later
#       ser   smp   tr      value
#  1: ser 1 smp 1 tr 1 0.53178276
#  2: ser 2 smp 1 tr 1 0.28574271
#  3: ser 3 smp 1 tr 1 0.62991366
#  4: ser 1 smp 2 tr 1 0.31073376
#  5: ser 2 smp 2 tr 1 0.36098971
# ---                            
# 20: ser 2 smp 1 tr 4 0.38049334
# 21: ser 3 smp 1 tr 4 0.14170226
# 22: ser 1 smp 2 tr 4 0.63719962
# 23: ser 2 smp 2 tr 4 0.07100314
# 24: ser 3 smp 2 tr 4 0.11864134

Observe cuán fácil fue esto, con todas nuestras etiquetas de dimensiones que se filtran en el formato largo. Una de las data.tables y data.tables de data.tables es la capacidad de realizar fusiones indexadas entre data.table s (al igual que las uniones indexadas de MySQL). Así que aquí, haremos eso para unir la class a nuestros datos:

dt <- dt.raw[J(trial.labs, class=trial.class)]  # on the fly mapping of trials to class
#          tr   ser   smp     value class
#  1: Trial 1 ser 1 smp 1 0.9648542     A
#  2: Trial 1 ser 2 smp 1 0.7285704     A
#  3: Trial 1 ser 3 smp 1 0.3142587     A
#  4: Trial 1 ser 1 smp 2 0.4134501     A
#  5: Trial 1 ser 2 smp 2 0.1393077     A
# ---                                    
# 20: Trial 4 ser 2 smp 1 0.2432201     A
# 21: Trial 4 ser 3 smp 1 0.2671743     A
# 22: Trial 4 ser 1 smp 2 0.5160964     A
# 23: Trial 4 ser 2 smp 2 0.7702306     A
# 24: Trial 4 ser 3 smp 2 0.8568685     A

Algunas cosas para entender:

1. J crea una data.table de data.table de vectores
2. intentando subconjuntar las filas de una tabla de datos con otra tabla de datos (es decir, utilizando una tabla de datos como el primer argumento después de la llave en [.data.table ) hace que data.table a la izquierda (en el lenguaje de MySQL) la tabla externa ( dt en este caso) a la tabla interna (la creada sobre la marcha por J ) en este caso. La unión se realiza en la (s) columna (s) de la key de la data.table externa, que como habrán notado, la definimos en el paso de conversión melt / data.table anterior.

Deberá leer la documentación para comprender completamente lo que está sucediendo, pero piense que J(trial.labs, class=trial.class) es efectivamente equivalente a crear una data.table con data.table(trial.labs, class=trial.class) , excepto que J solo funciona cuando se usa dentro de [.data.table .

Entonces, en un solo paso, tenemos nuestros datos de clase unidos a los valores. Nuevamente, si necesita álgebra matricial, opere primero en su matriz y luego, en dos o tres comandos sencillos, vuelva al formato largo. Como se señaló en los comentarios, es probable que no desee ir y venir del formato largo al de matriz, a menos que tenga una buena razón para hacerlo.

Una vez que las cosas están en data.table , puede agrupar / agregar sus datos (similar al concepto de estilo split-apply-combine) con bastante facilidad. Supongamos que queremos obtener estadísticas de resumen para cada class : combinación de sample :

dt[, as.list(summary(value)), by=list(class, smp)]

#    class   smp    Min. 1st Qu. Median   Mean 3rd Qu.   Max.
# 1:     A smp 1 0.08324  0.2537 0.3143 0.4708  0.7286 0.9649
# 2:     A smp 2 0.10130  0.1609 0.5161 0.4749  0.6894 0.8569
# 3:     B smp 1 0.14050  0.3089 0.4773 0.5049  0.6872 0.8970
# 4:     B smp 2 0.08294  0.1196 0.1562 0.3818  0.5313 0.9063

Aquí, simplemente le damos a data.table una expresión ( as.list(summary(value)) ) para evaluar para cada class , smp subconjunto de los datos (como se especifica en la expresión by ). Necesitamos as.list para que los resultados sean data.table por data.table como columnas.

También podría haber calculado fácilmente momentos (por ejemplo, list(mean(value), var(value), (value - mean(value))^3 ) para cualquier combinación de las variables clase / muestra / prueba / serie.

Si quieres hacer transformaciones simples a los datos, es muy fácil con data.table :

dt[, value:=value * 10]  # modify in place with `:=`, very efficient
dt[1:2]                  # see, `value` now 10x    
#         tr   ser   smp    value class
# 1: Trial 1 ser 1 smp 1 9.648542     A
# 2: Trial 1 ser 2 smp 1 7.285704     A

Esta es una transformación in situ, por lo que no hay copias de memoria, lo que hace que sea rápido. Generalmente data.table intenta usar la memoria de la manera más eficiente posible y, como tal, es una de las formas más rápidas de realizar este tipo de análisis.

Trazando desde un formato largo

ggplot es fantástico para trazar datos en formato largo. No entraré en los detalles de lo que está sucediendo, pero espero que las imágenes te den una idea de lo que puedes hacer

library(ggplot2)
ggplot(data=dt, aes(x=ser, y=smp, color=class, size=value)) + 
  geom_point() +
  facet_wrap( ~ tr)

ggplot(data=dt, aes(x=tr, y=value, fill=class)) + 
  geom_bar(stat="identity") +
  facet_grid(smp ~ ser)

ggplot(data=dt, aes(x=tr, y=paste(ser, smp))) + 
  geom_tile(aes(fill=value)) + 
  geom_point(aes(shape=class), size=5) + 
  scale_fill_gradient2(low="yellow", high="blue", midpoint=median(dt$value))

Tabla de datos -> Matriz -> Tabla de datos

Primero, debemos volver a acast (desde el paquete reshape2 ) nuestra tabla de datos a una matriz:

arr.2 <- acast(dt, ser ~ smp ~ tr, value.var="value")
dimnames(arr.2) <- dimnames(arr)  # unfortunately `acast` doesn''t preserve dimnames properly
# , , tr = Trial 1
#        smp
# ser        smp 1    smp 2
#   ser 1 9.648542 4.134501
#   ser 2 7.285704 1.393077
#   ser 3 3.142587 1.012979
# ... omitted 3 trials ...

En este punto, arr.2 ve igual que arr , excepto con valores multiplicados por 10. Tenga en cuenta que tuvimos que abandonar la columna de la class . Ahora, hagamos algo de álgebra matricial trivial

shuff.mat <- matrix(c(0, 1, 1, 0), nrow=2) # re-order columns
for(i in 1:dim(arr.2)[3]) arr.2[, , i] <- arr.2[, , i] %*% shuff.mat

Ahora regresemos al formato largo con melt . Tenga en cuenta el argumento key :

dt.2 <- data.table(melt(arr.2, value.name="new.value"), key=c("tr", "ser", "smp"))

Finalmente, volvamos a unir dt y dt.2 . Aquí debes tener cuidado. El comportamiento de data.table es que la tabla interna se unirá a la tabla externa en función de todas las claves de la tabla interna si la tabla externa no tiene claves. Si la tabla interna tiene claves, data.table se unirá a la clave. Este es un problema aquí porque nuestra tabla externa prevista, dt ya tiene una clave solo tr desde anterior, por lo que nuestra unión solo se realizará en esa columna. Por eso, debemos soltar la tecla en la tabla externa o restablecer la clave (elegimos la última aquí):

setkey(dt, tr, ser, smp)
dt[dt.2]
#          tr   ser   smp    value class new.value
#  1: Trial 1 ser 1 smp 1 9.648542     A  4.134501
#  2: Trial 1 ser 1 smp 2 4.134501     A  9.648542
#  3: Trial 1 ser 2 smp 1 7.285704     A  1.393077
#  4: Trial 1 ser 2 smp 2 1.393077     A  7.285704
#  5: Trial 1 ser 3 smp 1 3.142587     A  1.012979
# ---                                             
# 20: Trial 4 ser 1 smp 2 5.160964     A  5.867905
# 21: Trial 4 ser 2 smp 1 2.432201     A  7.702306
# 22: Trial 4 ser 2 smp 2 7.702306     A  2.432201
# 23: Trial 4 ser 3 smp 1 2.671743     A  8.568685
# 24: Trial 4 ser 3 smp 2 8.568685     A  2.671743

Tenga en cuenta que data.table realiza uniones al hacer coincidir columnas clave , es decir, al hacer coincidir la primera columna de clave de la tabla externa con la primera columna / clave de la tabla interna, la segunda con la segunda, y así sucesivamente, sin considerar la columna nombres (hay un FR here ). Si sus tablas / claves no están en el mismo orden (como fue el caso aquí, si lo notó), o necesita volver a ordenar sus columnas o asegurarse de que ambas tablas tengan claves en las columnas que desea en el mismo orden ( lo que hicimos aquí). La razón por la que las columnas no estaban en el orden correcto es debido a la primera unión que hicimos para agregar la clase, la cual se unió en tr e hizo que esa columna se convirtiera en la primera en la data.table .