¿Cómo distanciar y trazar en un gráfico MDS los objetos según su forma compleja?

Question

Supongamos que tengo cuatro formas básicas de señal (azul, morado, rojo, verde). También he creado formas de transición entre ellas. Si observas detenidamente la imagen de abajo, puedes ver que por ejemplo la señal azul (A1) se transforma lentamente en morado (A5) - horizontalmente, o en rojo (E1) - verticalmente. En general, cuanto más cercana sea la diferencia entre dos señales, más similares son.

**

Estoy buscando algún método/algoritmo/técnica que sea capaz de:

1. Extraer la mayor cantidad de información posible sobre la complejidad de las señales
2. Mapear las señales en 2D (o 3D) de acuerdo a su similitud

Enlace a los datos fuente: cada señal está codificada en una imagen JPG (25 archivos jpg, cada uno de 100x180 píxeles) https://www.dropbox.com/sh/ynidcsjdymrh85f/AAAHVtYSG0GUvX3CEdDWWG42a

---

He pasado algún tiempo intentando resolver este problema, así que aquí agrego mi enfoque, que no produce el resultado deseado:

En primer lugar, he establecido el directorio de trabajo, he cargado todos los 25 archivos jpg a una lista de R...

setwd("añadir directorio donde descargaste los archivos jpg aquí")

# librerías requeridas
library(jpeg)
library(raster)
library(asbio)

sgnl.vctr<-c("A1.jpg","A2.jpg","A3.jpg","A4.jpg","A5.jpg",
             "B1.jpg","B2.jpg","B3.jpg","B4.jpg","B5.jpg",
             "C1.jpg","C2.jpg","C3.jpg","C4.jpg","C5.jpg",
             "D1.jpg","D2.jpg","D3.jpg","D4.jpg","D5.jpg",
             "E1.jpg","E2.jpg","E3.jpg","E4.jpg","E5.jpg") 

sgnl.list <- list()      # sgnl.list contiene todas las 25 señales
for (i in 1:length(sgnl.vctr)){
sgnl.list[[i]] <- readJPEG(sgnl.vctr[i])
}                                                       

Experimenté problemas con los valores de píxeles (rango de 0 a 1), por lo tanto los reclasifiqué en binarios (0 o 1).

# reclasificación de valores
for (i in 1:25) {
sgnl.list[[i]][1:100,1:180,1][sgnl.list[[i]][1:100,1:180,1] > 0.5]<- 2
sgnl.list[[i]][1:100,1:180,1][sgnl.list[[i]][1:100,1:180,1] <= 0.5]<- 1
sgnl.list[[i]][1:100,1:180,1][sgnl.list[[i]][1:100,1:180,1] == 2]<- 0
}

Luego, extraje vectores binarios (0-1) de cada archivo jpeg. Si alguien sabe cómo acortar el procedimiento abajo, por favor edite el código R.

# Fila A
A1<-as.vector(sgnl.list[[1]][1:100,1:180,1])
A2<-as.vector(sgnl.list[[2]][1:100,1:180,1])
A3<-as.vector(sgnl.list[[3]][1:100,1:180,1])
A4<-as.vector(sgnl.list[[4]][1:100,1:180,1])
A5<-as.vector(sgnl.list[[5]][1:100,1:180,1])
# Fila B
B1<-as.vector(sgnl.list[[6]][1:100,1:180,1])
B2<-as.vector(sgnl.list[[7]][1:100,1:180,1])
B3<-as.vector(sgnl.list[[8]][1:100,1:180,1])
B4<-as.vector(sgnl.list[[9]][1:100,1:180,1])
B5<-as.vector(sgnl.list[[10]][1:100,1:180,1])
# Fila C
C1<-as.vector(sgnl.list[[11]][1:100,1:180,1])
C2<-as.vector(sgnl.list[[12]][1:100,1:180,1])
C3<-as.vector(sgnl.list[[13]][1:100,1:180,1])
C4<-as.vector(sgnl.list[[14]][1:100,1:180,1])
C5<-as.vector(sgnl.list[[15]][1:100,1:180,1])
# Fila D
D1<-as.vector(sgnl.list[[16]][1:100,1:180,1])
D2<-as.vector(sgnl.list[[17]][1:100,1:180,1])
D3<-as.vector(sgnl.list[[18]][1:100,1:180,1])
D4<-as.vector(sgnl.list[[19]][1:100,1:180,1])
D5<-as.vector(sgnl.list[[20]][1:100,1:180,1])
# Fila E
E1<-as.vector(sgnl.list[[21]][1:100,1:180,1])
E2<-as.vector(sgnl.list[[22]][1:100,1:180,1])
E3<-as.vector(sgnl.list[[23]][1:100,1:180,1])
E4<-as.vector(sgnl.list[[24]][1:100,1:180,1])
E5<-as.vector(sgnl.list[[25]][1:100,1:180,1])

Los vectores fueron comparados con la función Kappa para obtener el valor total de acuerdo. Ver este enlace: https://stackoverflow.com/questions/24534192/how-to-compare-all-possible-combinations-of-objects-in-r-by-loop/24534794#comment37992299_24534794 (Muchas gracias a @digEmAll)

# bucle de:
#            Kappa(sgnl.list[[1]][1:100,1:180,1],
#                  sgnl.list[[2]][1:100,1:180,1])$ttl_agreement
M <- rbind(A1,A2,A3,A4,A5,
           B1,B2,B3,B4,B5,
           C1,C2,C3,C4,C5,
           D1,D2,D3,D4,D5,
           E1,E2,E3,E4,E5)

res <- outer(1:nrow(M),
             1:nrow(M),
             FUN=function(i,j){
               # i y j son 2 vectores de la misma longitud que contienen 
               # las combinaciones de los índices de las filas. 
               # por ejemplo (i[1] = 1, j[1] = 1) (i[2] = 1, j[2] = 2)) etc...
               sapply(1:length(i),
                      FUN=function(x) Kappa(M[i[x],],M[j[x],])$ttl_agreement )
             })

row.names(res) <- c("A1","A2","A3","A4","A5",
                    "B1","B2","B3","B4","B5",
                    "C1","C2","C3","C4","C5",
                    "D1","D2","D3","D4","D5",
                    "E1","E2","E3","E4","E5")
colnames(res) <- c("A1","A2","A3","A4","A5",
                   "B1","B2","B3","B4","B5",
                   "C1","C2","C3","C4","C5",
                   "D1","D2","D3","D4","D5",
                   "E1","E2","E3","E4","E5")

Finalmente, la matriz de similitud (objeto res) se utilizó para escalamiento multidimensional...

# mds basado en la matriz ttl_agreement
d <- as.dist(res)
mds.coor <- cmdscale(d)
plot(mds.coor[,1], mds.coor[,2], type="n", xlab="", ylab="")
text(jitter(mds.coor[,1]), jitter(mds.coor[,2]),
     rownames(mds.coor), cex=0.8)
abline(h=0,v=0,col="gray75")

Sin embargo, como puedes ver (gráfico de la izquierda), las cuatro señales básicas no se separaron como esperaba. ¿Alguien conoce una mejor solución que conduzca al resultado deseado (gráfico de la derecha)?

**

Answer 1

Esta puede ser solo una respuesta parcial porque no creo que la trama que esperas sea realmente lo que está en los datos, especialmente la "paralelidad y continuidad" de las señales intermedias. Especularé sobre las razones de eso a continuación.

Pero creo que pude llegar a lo que estás buscando en términos de las cuatro señales basales A1, A5, E1, E5. Es decir que yacen en el borde del espacio embebido, que las señales opuestas están más o menos diametralmente opuestas (es decir, A1, E5 y E1, A5 respectivamente) y que las señales vecinas se conservan (así A1, A5 y E1, E5 respectivamente).

En general, creo que el MDS estándar (es decir, con una matriz de pesos de entrada compuesta solo por 1) realmente no te da lo que quieres porque en realidad estás buscando una reducción de dimensionalidad no lineal que tenga alguna característica de localización, es decir, que las "distancias" más grandes no deberían ser preservadas, sino las distancias locales deberían serlo. Hay una variedad de algoritmos que hacen eso. Uno bastante popular que se usa a menudo para casos como el tuyo se llama t-SNE o Embedding Estocástico de Vecinos Distribuidos en t. En la página enlazada encontrarás bastante información.

Calculé el embedding t-SNE para tu distancia kappa. Ten en cuenta que uso una alta perplejidad para forzar una forma que se parezca a la que buscas. Se puede obtener por

set.seed(1)
library(tsne)
tsne.coor1 <- tsne(res,perplexity=25)
rownames(tsne.coor1) <- c("A1","A2","A3","A4","A5",
                   "B1","B2","B3","B4","B5",
                   "C1","C2","C3","C4","C5",
                   "D1","D2","D3","D4","D5",
                   "E1","E2","E3","E4","E5")
plot(tsne.coor1[,1], tsne.coor1[,2], type="n", xlab="", ylab="")
text(tsne.coor1[,1], tsne.coor1[,2],
     labels=row.names(tsne.coor1), cex=0.8)
abline(h=0,v=0,col="gray75")

y se ve así

Las señales básicas están etiquetadas en rojo. Como puedes ver, la estructura de similitud capturada por tu preprocesamiento y medida de distancia no sugiere que las señales intermedias sean realmente como esperas. Por ejemplo, el camino de A1 a E1 tiene señales intermedias E4, D5, A3 y no B1 a través de D1. ¡Pero esto es lo que tus datos te dicen! Entonces, la similitud entre las señales que están dentro del cascarón convexo del espacio incrustado sugiere que el patrón claro no se conserva.

Hay dos explicaciones obvias:

1. Se pierde algo de información al mapear a bajas dimensiones. La información que se pierde podría ser la que estabas buscando realmente.
2. La medida de distancia podría no captar lo que realmente te importa. Estoy de acuerdo con @ttnphns en esto.

[\Editar] (¡gracias @ttnphns!) Para investigar el último punto, probé otras medidas de similitud para matrices binarias. Para esta alta perplejidad no llevó a resultados discerniblemente diferentes en la forma, pero sí lo hizo en los arreglos (usé un modelo de similitud de gravitación y una similitud binaria asimétrica, la que devuelve dist(x,method="binary")). Para la baja perplejidad, los efectos en la visualización con la distancia kappa son pequeños, pero para las otras distancias no lo son. Para ilustrar los resultados con distancias asimétricas:

res <- dist(M,method="binary")
set.seed(1)
tsne.coor2 <- tsne(res,perplexity=25)
tsne.coor3 <- tsne(res,perplexity=3)
rownames(tsne.coor2)<-rownames(tsne.coor3) <- c("A1","A2","A3","A4","A5",
                   "B1","B2","B3","B4","B5",
                   "C1","C2","C3","C4","C5",
                   "D1","D2","D3","D4","D5",
                   "E1","E2","E3","E4","E5")
par(mfrow=c(1,2))
plot(tsne.coor2[,1], tsne.coor2[,2], type="n", xlab="", ylab="")
text(tsne.coor2[,1], tsne.coor2[,2],
     labels=row.names(tsne.coor3), cex=0.8)
abline(h=0,v=0,col="gray75")
plot(tsne.coor3[,1], tsne.coor3[,2], type="n", xlab="", ylab="")
text(tsne.coor3[,1], tsne.coor3[,2],
     labels=row.names(tsne.coor2), cex=0.8)
abline(h=0,v=0,col="gray75")

y aquí están los resultados

Por lo tanto, para alta perplejidad, la diferencia en la forma de la proyección de las señales es bastante similar para ambas medidas de distancia. Para baja perplejidad, los resultados cambian. Ten en cuenta que aunque se ve diferente de lo que originalmente pretendías al usar la distancia binaria asimétrica, la transición desde la señal base a través de los otros estados parece estar mejor preservada en la gráfica de baja perplejidad, especialmente en columnas! . Esto hace que la especulación 2 de usar la distancia kappa como una medida de distancia no muy adecuada sea más probable. Puedes probar las distancias dist(x, method="binary") o cluster::daisy(x, metric= "gower")) (que es el coeficiente de dice, creo).

[\Fin Edit]

Ten en cuenta que t-SNE tiene inicializaciones aleatorias, por lo que puede verse un poco diferente en tu caso --- no estoy seguro de si set.seed tiene un efecto en la implementación de R. La inicialización aleatoria es algo que podría acercarte más a lo que necesitas de todos modos. Como lo explican los autores:

En contraste con, por ejemplo, PCA, t-SNE tiene una función objetivo no convexa. La función objetivo se minimiza utilizando una optimización de descenso de gradiente que se inicia de manera aleatoria. Como resultado, es posible que diferentes ejecuciones te den soluciones diferentes. Ten en cuenta que está perfectamente bien ejecutar t-SNE varias veces (con los mismos datos y parámetros) y seleccionar la visualización con el valor más bajo de la función objetivo como tu visualización final.

Puedes jugar con las distancias mencionadas anteriormente y algunos de los parámetros, especialmente la perplejidad, para quizás acercarte más a lo que necesitas en una dirección u otra. Espero que esto sea de ayuda como un primer paso.