11 respuestas en este tema

[Delphius]

Hola,
Como dice el título, ¿alguien me puede dar una manito en interpretar un poco de código escrito en Matlab?

Verán tengo a modo de ejemplo cierto algoritmo para adaptar a Delphi y hay algo que no me cuadra o me cuesta entenderle. De todo el gran algoritmo lo que me cuesta entender es como juega el tema del operador : (que yo vulgarmente llamo iterador) cuando uno lo combina con el operador de traspuesta (') y en una multiplicación de matrices.
Tengo una matriz u de tamaño N x M, una matriz A de NxM previamente reservada. Y lo que sigue es aplicar una operación que me cuesta entenderle, porque todo lo anterior ya lo tengo muy bien armado.

El lio está en esto:

delphi

w=zeros(M,M);
for (i=1:M)
   for (j=1:M)
      w(i,j)=u(:,j)'*A(:,j);
   end
end

Justamente lo que me cuesta entenderle es la operación real que sucede en el cálculo de w(i,j). En teoría se que esto representa una operación punto a punto entre los elementos de u y A. Pero... el que se aplica una traspuesta me induce a pensar si lo que el u original es una matriz NxM ahora será MxN. De se así entonces si el iterador : está en la fila ¿al aplicar transpuesta será sobre la columna? ¿O es que itera sobre la fila de la matriz traspuesta?
¿Que tiene prioridad?

Si alguien le entiende esa partecita se lo agradecería. Porque cada vez que lo veo me mareo. Quizá sea una tontera pero tengo la cabeza ya quemada entre tantas cosas que he visto en el día.

EDITO:
Por si sirve de ayuda, se que la matriz resultante debe ser MxM. La parte teórica de mis apuntes dice justamente esto:

w_k = u_k^T*A, for k = 1 to 1,..., M

El texto que acompaña a esto: "Esto describe un conjunto de multiplicaciones y sumatorias punto-a-punto. Los pesos w forman un vector Omega^T=[w₁, w₂, ..., w_M]"

Inicialmente ya tenía código escrito en Delphi que en esta última parte no me cuadraba y tuve que recurir a material complementario que se da en Matlab. Y rompió mis moldes, y complicó el escenario final en como se debe leer cada vector Omega_k.

Saludos,

[Delphius]

Mi primer análisis, aunque sospecho que no es del todo correcto, me llevó a lo siguiente:
Asumiendo que primero efectúa una traspuesta, la matriz u de NxM pasa a ser MxN por lo que una referencia u(:,j) pasa a ser en un escenario traspuesto una matriz hipotética uT(j,.
Luego ya que en realidad este operador : implica iterar por esta dimensión y como ésta siempre es de orden N, entonces u(n,m) equivale a uT(m,n).

Como La traspuesta sólo implica cambiar filas por columnas (o viceversa) y los valores no se alteran, la posición (i,j) en u para su traspuesta equivale a (j,i).

Además este operador : debe coincidir en dimensión con el de A. Y efectivamente conciden, ya que A es de NxM.
Entonces me dije, bueno, tal operacíon producto interno entre u y A va a ser algo como MxN * NXM. Y lo que haría es iterar sobre N.

En resumen mi primer propuesta fue algo como:

delphi

for i := 1 to M do
  for j := 1 to M do
    for k = 1 to N do
      w(i,j) := w(i,j) + u(k,j) * A(k,i)

El asunto es que cuando reviso mi literatura, bien señala que para un vector columna de A (Nx1), un A_i llamemosle, en realidad se tiene esta operación:

w_k = u_k^TA_i, para k = 1...M

De modo que estos m valores, w₁...w_m se deben agrupar en un vector Omega_i^T=[w₁, ... w_m]. Como este vector está traspuesto en realidad se trata de un vector columna.
Luego se debe repetir este mismo procedimiento para todos los vectores columnas i = 1...M de A.

Esto da origen a esa matriz w(i,j) que intento calcular. Con el algoritmo que he propuesto en realidad llego a un escenario traspuesto. Aparentemente en lugar de obtener un vector omega columna, obtendría uno fila.

Pensé entonces que me bastaría con invertir los índices de w() en el código. Pero he aquí que el código de muestra que acompaña a la teoría que acabo de dar expresa justamente W(i,j) y no W(j,i). Y Me extrañaría un error en el código cuando en el resto de todo el documento las cosas están perfectamente bien.
He extrapolado el 90% del código y efectivamente funciona. En ese 90% hay operaciones : pero no se ha requerido de estar aplicando una traspuesta en combinación de éste.

Y eso ya me hizo extraño.

¿Algún norte de como interpretar esto?

¿Donde está Sergio cuando lo necesito? 

Saludos,

[Sergio]

Aqui llega Sergio!

Casi lo tienes, estas multiplicando 2 matrices, solo eso, pero en matlab la forma de escribirlo es al minimo esfuerzo, y cuesta desenmarañarlo.

Si tienes dos matrices, U de NxM y A de MxN, al multiplicar UxA obtienes una matriz W que tendra tantas filas como la primera matriz (N filas tiene U) y tantas columnas como columnas tenga la segunda (A tiene N columnas), asi que en este caso, W es una matriz cuadrada de NxN.

Cada elemento de W, w(i,j), es el resultado de multiplicar la fila i de U por la fila j de A, y para poderse multiplicar esas filas/columnas, han de tener el mismo numero de elementos: M en este caso.

Asi que:

w(i,j)= Fila i de U * Columna j de A = Sumatorio para k de 1 a M de U(i,k) * A(k,j)

Es decir, que casi casi lo tenias ya:

delphi

for i:= 1 to N do begin
  for j:= 1 to N do begin
    w(i,j):= 0;
    for k:= 1 to M do begin
      w(i,j):= w(i,j) + U(i,k) * A(k,j);
    end;
  end;
end;

[Delphius]

Gracias Sergio por venir al rescate, aunque en lugar de encontrarte con una bella dama va a ser un tipo de 30, gordo y medio pelado. 

En realidad A es de NxM. De las mismas dimensiones que U.
Por lo que la teoría dice que tal multiplicación no es posible: (NxM) x (NxM)

El código Matlab, acompañando en cierta forma la teoría tiene en su escritura justamente esto:

W(i,j) = U(:,j)'*A(:,i)

De lo que yo entiendo es que para que tenga sentido tal multiplicación de producto interno es que realiza una traspuesta sobre U, de modo que lo que era NxM sea MxN. Y entonces lo que era U(:,j) ahora sea un UT(j, Y allí si ya entonces puedo iterar por k=1..N Tal como dices.

A mi lo que me confunde es justamente como entenderle esa vuelta de tuerca.

Obviando un poco el codigo matlab, los textos aseguran que lo que debe hacerse es esto:

w_k = u_k^TA_i, para k=1..M

El texto que acompaña esta formula dice: This describes a set of point-by-point image multiplications and sumations.

Siendo A_i justamente un vector columna de A: Nx1.
Estos w₁..w_m se organizan en un vector O^T = [w₁..w_m]. Por lo que el O original es una columna Mx1. Al generar cada columna O_i se consigue esa matriz W de MxM.

Para que veas de donde sale U, inicialmente, viene de esto:

delphi

zeros(N:M)
for (i=1:M)
  for (j=1:M)
    U(:,i) = U(:,i) + V(i,j) * A(:,j)

V es una matriz de MxM autovectores calculados desde una matriz de covarianza L (MxM) = A^TA. A es el conjunto de datos centrados, con media cero. Es decir que A = A - Media.

Ya con la matriz W y U devo "volverme sobre mis pasos", y obtener la reconstrucción de A.
Areconstruida_i = sumatoria (W_k*U_k) siendo w₁..w_m ese vector O_i.

Seguramente tu te haces una idea del algoritmo que estoy aplicando. 

Por eso me es importante lograr determinar como es que realmente está siendo calculado W. Y a estas horas de mal dormir y de pensar ya dudo de mi y del código Matlab. No se si pensar que el código está bien y yo debo tomar justamente tal como dice la teoria Y hacer que cada vector columna de W sea el vector O. O si es que desafiando a la propia teoría el código está generando la matriz W e interpretar a cada fila de W como el vector O.

¿Porqué digo esto? Por esto:
The weights form a vector O^T = [w₁..w_m] that describes the contribution of each eigenface (NOTA: Eigenface es cada columna U_k) in representing the input image, treating the eigenfaces as basis set for face images.

Y ahora cuando veo tu código me quedé en modo pensativo.

Espero no haberte confundido a vos ahora. 

Saludos,

[Sergio]

Es que si quieres hablar mas del significado de las matrices que de la operacion... eso se avisa!

Traduciendo lo que dices (si es que lo entiendo) andas detectando caras a partir de fotos.

Tienes M caras almacenadas con N propiedades extraidas de ellas (distancia ojos etc), y cuando le metes tu cara nueva, obtienes los w que son pesos, es decir, te pareces a la cara 1 en un w1 por ciento, etc.

Y ahora tienes los w pero necesitas la matriz A con los valores de esa cara concreta (Aij seria: como de parecido eres con la cara Mi en la propiedad j).

Explicame esto, lo que haces de verdad y lo que es cada cosa, las formulas son secundarias.

Por cierto, si no comentas que la cosa va de "caras" no lo pillo nunca, ser matematico no te hace adivino y nunca he trabajado en reconocimiento facial.

Sobre si hacer caso a teoria o ejemplo que funciona: Si el ejemplo funciona esta ajustado a la teoria, copia el codigo y luego te piensas porque esta ajustado y no lo ves.

"En teoria, teoria y la practica son lo mismo. En la practica, no." 

[Delphius]

Asi es amigo. Es una aplicación práctica del método estadístico PCA, O Análisis de Componentes Principales con reducción de dimensión.

Básicamente es armar una matriz A de NxM siendo N el tamaño de las imagenes (N = ancho x alto). y M la cantidad de fotos. A este conjunto se le calcula un vector columna Nx1 la media y se procede a restar a cada columna de esta matriz A el vector Media:
A = A - Media

Luego hay que proceder a armar una matriz de covarianza. Si hiciéramos la matriz original C = AA^T Dicha matriz C sería de NxN. Lo que es algo impracticable. En su lugar se aplica reducción de dimensión. Esto se consigue haciendo lo opuesto: A^TA, obteniendo una matriz de covarianza L de MxM.
A esta matriz L se calcula los autovectores y autovalores, que se almacenan en V (MxM) y Lambda (1xM) respectivamente.

El "truco" está en que ahora estos autovectores deben ser llevados a la dimensión original. Por teorema de álgebra esto se obtiene así:

U = sumatoria(V_ij*A_j)

Siendo U una matriz de NxM. Donde Cada vector columna U_j j = 1..M será una autocara o eigenface. Una autocara es una combinación lineal de todas las caras. Para volver a recrear las caras hace falta la matriz de pesos, que determina los pesos o ponderación de cuanto contribuye cada autocara para obtener la original. De esta forma se tendrá tantos vectores O_i = w₁..w_m como caras que existan.

La matriz W tendrá una dimensión MxM. Tal fórmula ya la expresé, y es la que en Matlab en "cierta forma" simula.
Todo mi desarrollo anterior, y posterior a lo que sigue a W coincide tanto en teoría y práctica. Más es el armado de la matriz W en donde aparentemente hay cierto dezliz, o bien en la teoría o en la práctica... o en ambos.

Cuando ya se tiene la matriz de pesos W, y dada una nueva imagen llamemosla F (Nx1) a la que procedemos a restarle la media (F = F - Media) procedemos a obtener sus propios pesos al projectarla sobre las autocaras:

w_kf = U_k^TF

Con los pesos w_kf agrupados en un O_f ahora vemos que tan distante está respecto a cada O_i. Esto no es más que calcular la distancia d(Of, Oi). Aquel O_i que ofrezca la menor distancia y cumpla con cierto valor de umbral de aceptación correspondería al identificado.

Luego a fin de comprobar esto generamos una reconstrucción de esa i-ésima cara:
R_i = sumatoria(w_k * U_k) siendo los pesos wk los correspondientes al O_i. De aquí es que es fundamental organizar cada w(i,j) de manera apropiada.

Si resulta ser que d(F, Ri) también cumple con el valor de umbral efectivamente es un reconocimiento positivo. Y ya de yapa para reconstruir la imagen original sin central basta con hacer: R_i + Media

Volviendo al tema del W, lo que cuesta vislumbrar es que procedimiento aplicar. Y con más razón entenderle a ese pedacito de código Matlab. Puedes ver el paper si deseas aqui por si te interesa y te ayuda a comprender mejor.

Otras fuentes que he consultado se pasan justamente al armado de W a sus anchas... algunos aseguran que O es columna, otros lo ven como fila. Aún cuando citan a este mismísimo documento. Unos pocos exponen que para el armado de w no se aplica traspuesta  e incluso éstos también divagan en si W debe leerse por columna o fila. Es de locos. 

Saludos,

[Delphius]

También dejo el paper que toma el trabajo original y expone el respectivo código matlab. Puede ser consultado aquí.

Saludos,

[Delphius]

Mi cabeza hace esfuerzo de lógica, y me inclina hacia que hay una pifia en la muestra de Matlab. Y debe revisarse dicho código para ser congruente con la formación de la matriz de pesos.

Como el producto interno entre dos vectores columnas A y B se define como: A^TB. Y en última todo se resume en una simple sumatoria de los productos de sus componentes: A₁*B₁ + A₂*B₂ + ... No hace falta otra cosa especial, ni interesa realmente si están en columna o fila.

Mi actual propuesta es la siguiente:

delphi

for i := 1 to M do
  for j := 1 to M do
    for k := 1 to N do
      W(j, i) := W(j, i) + U(k, j) * A(k, i);

De esta forma consigo dos cosas:
1) No tener que malgastar en una nueva matriz traspuesta cuando se puede obtener los datos de la propia
2) Expresar los m pesos w_k asociados al vector columna A_i como sugiere tanto la teoría en un vector O_i columna de Mx1. Y de esta forma el vector O_i obtenido respondería a la contribución de cada eigenface respecto a la imagen i-ésima.

Si siguiera mi propuesta inicial, tal como está en código Matlab la matriz W estaría traspuesta. Y de seguir la propia teoría, si se agrupara los pesos en forma de columna se estaría leyendo cualquier cosa (Vendría a ser algo como la "contribución" del eigenface i-ésimo a todas las imágenes). Para arreglar las cosas en lugar de hacer Un O_i columna debiera ser fila.

Quedo a la espera de Sergio a ver el dictamen final. Y queda como incógnita esa supuesta pifia en la muestra de matlab.

Saludos,

[Sergio]

He tenido que leerlo dos veces para seguirte... tendre que meterme en esto algun dia, suena endiabladamente complicado pero elegante!

Creo que te estas liando en ese detalle de como armar w: solo tienes una manera de hacerlo, fijate en los for, no puedes cambiar nada alli, y con esos for, solo puedes multiplicar esos elementos de cada amtriz, no puedes intercambiar los indices ni nada por las dimensiones de las matrices... solo se puede multiplicar las dos matrices!

Si comparas tu propuesta final con la mia, es la misma exactamente: estas multiplicando 2 matrices, y eso solo se puede hacer de una manera.

Que w lo veas como fila o columna, da igual, solo te tocara trasponerla o no, detalles, si el paper dice fila y mathlab columna, da igual, es lo mismo.

Por cierto, en tu codigo te falta inicializar los W a cero como en el mio, no te olvides.

De todas formas, pasa de ese detalle y usa multiplicacion de matrice,s termina tu codigo, y hazle unas pruebas: si funciona, w esta bien armado. No tendras otra forma de convencerte que probandolo (yo hago eso cuando la formulacion se me complica ya demasiado).

[Delphius]

No se puede hacer la multiplicación de matrices de toda la vida.  Tanto U como A tienen dimensión NxM. Por tanto: (NxM) x (NxM) no es posible.
De la única forma en que es posible tal cosa es aplicando traspuesta a U, o bien a A. Aplicando a U queda: UT: (MxN). Y allí si vale:

UT*A = (MxN) x (NxM) = MxM)

Tu te confundiste en tu primer post.

Luego tengo:
Indice i: 1..M >> I-ésimo espacio de cara
Indice j: 1..M >> J-ésima autocara
Indice k: 1..N >> Itera sobre los elementos interno de la autocara y la cara. Equivale a iterar por los pixeles

Para cada espacio de cara se debe obtener los M pesos asociados de las M autocaras.

Ya que leer una UT equivale a leer la misma U con índices invertidos, me parece válido hacer eso. Es decir que para cualquier valor j=1..M y k=1..N UT(j,k) = U(k, j).

Internamente cada posición de matriz resultante C de una multiplicación A*B es un producto interno entre la fila A y la columna B. Obviando de verlo como multiplicación de matrices de toda la vida y quedarme con la idea de un producto interno entre dos vectores, Nx1 y Nx1. El único requisito es que tengan el mismo "rango". Y entonces dije: bueno, me vale como estén ordenados realmente si todo puede resolverse iterando sobre sus índices adecuados.



Coincido en que al final pareciera ser que da lo mismo tener a W como filas o columnas. Después de todo otras referencia divagan en este punto. Lo importante es que al momento de leer cada vector O, se respete el orden correcto sino se leerá cualquier cosa.

La verdad es que en ese punto todos los documentos no son tan bien explicativos y lo dan de forma un tanto vaga. Al comienzo estaba como tu, viendo todo como si fueran meras multiplicaciones de matrices y buscaba la forma de que encajara todo... me pasó cuando intentaba hacer que U calzara como multiplicación de matrices hace un tiempo. Luego a medida que empecé a mirar de reojo el texto empiezo a notar que lo que aparentemente es algo no lo es. En el caso del cálculo de la matriz W quedaría mucho más claro si lo representaran como corresponde: w_ki = U_k^T·A_i

Viendo ese punto ya no te quedan dudas de que es un producto interno entre los vectores de las matrices.

No si es que Turk y cia, lo hicieron a drede para evitar que otros tomen sus ideas tan fácil o que. Lo que si es que te rompe el coco un buen rato.

Saludos,

[Sergio]

He visto cosas asi varias veces, multiplican matrices sin preocuparse de si antes hay que trasponer una de ellas, se da por hecho que lo haras porque de otra forma no sale nada, no se puede intentar, asi que ya dan por cierto que haras que los for no se "salgan de rango" poniendo (i,j) o (j,i) y se quedan tan panchos.

Es lo malo de que los algortimos los definan informaticos y fisicos, que no son matematicos! 

Pero bueno, al final funciona y eso es lo que importa.

[Delphius]

He visto cosas asi varias veces, multiplican matrices sin preocuparse de si antes hay que trasponer una de ellas, se da por hecho que lo haras porque de otra forma no sale nada, no se puede intentar, asi que ya dan por cierto que haras que los for no se "salgan de rango" poniendo (i,j) o (j,i) y se quedan tan panchos.

Es lo malo de que los algortimos los definan informaticos y fisicos, que no son matematicos! 

Pero bueno, al final funciona y eso es lo que importa.

¡Yo que me considero informatico matematiloco.. no te imaginas como me pongo!    Si leo una ecuación busco respetarla como corresponde. Se supone que para un trabajo de la altura de alguien que está en el MIT es lo mínimo que se pide.

Cuando veo que algo no cuadra, ya me agarra la locura. Esta parte me ha tenido mal. Esta madrugada me desperté soñando con eso, aunque ya de antes que vengo de un muy mal dormir. Mi cuerpo en general me pide descanso y asi no hay forma de concentrarse como corresponde.

Saludos,