2

Estoy construyendo una matriz (20x20) para unos indicadores que se tienen que emplear más adelante en un Balance Scorecard.

Para efectuar las estimaciones correspondientes, se toma la diagonal cuyos valores son unos rendimientos provenientes de históricos.

Con dicha diagonal se obtienen promedios para obtener una matriz triangular inferior: lower triangular matrix

Para que se visualice lo anterior aquí un pequeño esquema:

introducir la descripción de la imagen aquí

Dado que la matriz es relativamente pequeña (20x20) es válido emplear un ataque de fuerza bruta

promedioMatrizTriangular <- function(diagonal){
    len <- length(diagonal)
    m <- matrix(0, nrow=len, ncol = len)
    m <- diag(diagonal)

    for (j in 1:len) {        # columnas
        for (i in 1:len){     # filas

            if(i>j) m[i,j]<-mean(diagonal[j:i]) 

        }
    }

    round(m, digits = 1)
}

Corriendo la función anterior con la diagonal mostrada en el esquema previo, se tiene:

promedioMatrizTriangular(c(9,11,8,5))

     [,1] [,2] [,3] [,4]
[1,]  9.0  0.0  0.0    0
[2,] 10.0 11.0  0.0    0
[3,]  9.3  9.5  8.0    0
[4,]  8.2  8.0  6.5    5

Es decir funciona, pero... si se tuviera una matriz de 1000x1000 eso implicaría ejecutar if(i>j) un millón de veces.

Entonces para la respuesta, cualquier enfoque que se emplee y que conlleve a una reducción de ese millón de operaciones, será aceptada como valida.

Nota: Algunas aplicaciones prácticas de la matriz triangular:

3 respuestas 3

1

Consideraciones

En este caso la función promedioMatrizTriangular tiene un coste computacional de O(n²), esto no tiene un impacto significativo si la longitud de la diagonal (n) se mantiene en valores pequeños.

Al observarse la estructura de la matriz (4x4) se aprecia que no hay operaciones de mean en la última columna, eso implica que el indice j va de 1 a len-1.

Por otra parte, las filas sí tienen operaciones en la última de ellas, es decir: i=len.

Asimismo se proporciona un dato importante dentro del if: i>j esto sumado al hecho que los saltos en las filas de una matriz son discrecionales y como en este caso son de uno a uno, se puede afirmar que i=j+1.

Lo anterior se puede verificar de la siguiente manera:

Cuando j=1, el indice i comienza en la fila 2y cuando j=len-1, se tiene i=len.

Variante

Entonces por lo expuesto el código quedaría así:

promedioMatrizTriangular <- function(diagonal){
    len <- length(diagonal)
    m <- matrix(0, nrow=len, ncol = len)
    m <- diag(diagonal)

    for (j in 1:(len-1)) {        # columnas
        for (i in (j+1):len){     # filas

            m[i,j]<-mean(diagonal[j:i]) 

        }
    }

    round(m, digits = 1)
}

promedioMatrizTriangular(c(9,11,8,5))

     [,1] [,2] [,3] [,4]
[1,]  9.0  0.0  0.0    0
[2,] 10.0 11.0  0.0    0
[3,]  9.3  9.5  8.0    0
[4,]  8.2  8.0  6.5    5

Como sólo se recorre las celdas dónde se va ejecutar el mean, la cantidad de operaciones es de la forma: O(n(n-1)/2)

Resumiendo

Se puede graficar las dos O( ) -que se denotarán como f1 y f2- para ver su comportamiento cuando se tienen diagonales comprendidas en el rango 2:1000

f1 <- Vectorize(function(x) x*x, 'x')
f2 <- Vectorize(function(x) x*(x-1)/2, 'x')
r <- 2:1000

plot(r, f1(r), type = "l", lty=2,  col="red", yaxt = "n", ylab="Operaciones (Miles)", xlab="Longitud Diagonal")
lines(r, f2(r), type = "l",  col="blue")

# cambiando la nomenclatura del eje y: notese que en plot se limpio el eje y:  yaxt = "n"
axis(2, at = y <- seq(0, 1e6, by = 10e4),  labels = paste0(y/1000, "k"))
legend("topleft", col = c("Red","Blue"), lty=c(2,1), lwd=c(1,1), cex=0.75, c("Ataque Fuerza Bruta", "Indexando"))

introducir la descripción de la imagen aquí

1

Una forma de entender las matrices es verlas como un simple vector plano, de esa forma podemos estudiar mejor con que patrón se dan las combinaciones de los elemento de la diagonal a seleccionar, veamos:

introducir la descripción de la imagen aquí

Se puede ver que hay un patrón de como se combinan los números de la diagonal, en la parte de abajo esta indicada en cada caso, cuantos elementos de la diagonal hay que considerar y desde que posición de esta los tomaremos. El truco con esto, es generar los patrones para los vectores de cantidades y posiciones y asumir (me falta background matemático como para demostrarlo) que este patrón se repite con cualquier dimensión que elijamos (En mis pruebas lo hace).

promedioMatrizTriangular <- function(diagonal) {

    n = length(diagonal)
    # Genero los vectores para el recorte de la diagonal
    cantidades <- unlist(lapply(rep(n:1), FUN=function(x) c(rep(0,(n-x)),rep(1:x))))
    posiciones <- unlist(lapply(rep(1:n), FUN=function(x) c(rep(0,(x-1)), rep(x, n-(x-1)))))

    # Recorto la diagonal dado un índice 
    diag_mean <- function(x) {
        mean(diagonal[posiciones[x]:(posiciones[x] + (cantidades[x]-1))])
    }

    # Para cada valor de 1:n*n, recortamos la diagonal y calculamos la media
    # retornamos directamente la matriz
    matrix(lapply(1:(n*n), FUN=function(x) ifelse(posiciones[x]==0,0,diag_mean(x))), nrow=n)
}

Salida:

> promedioMatrizTriangular(c(9))
     [,1]
[1,] 9   
> promedioMatrizTriangular(c(9, 11))
     [,1] [,2]
[1,] 9    0   
[2,] 10   11  
> promedioMatrizTriangular(c(9, 11, 8))
     [,1]     [,2] [,3]
[1,] 9        0    0   
[2,] 10       11   0   
[3,] 9.333333 9.5  8   
> promedioMatrizTriangular(c(9, 11, 8, 5))
     [,1]     [,2] [,3] [,4]
[1,] 9        0    0    0   
[2,] 10       11   0    0   
[3,] 9.333333 9.5  8    0   
[4,] 8.25     8    6.5  5   
> promedioMatrizTriangular(c(9, 11, 8, 5, 6))
     [,1]     [,2] [,3]     [,4] [,5]
[1,] 9        0    0        0    0   
[2,] 10       11   0        0    0   
[3,] 9.333333 9.5  8        0    0   
[4,] 8.25     8    6.5      5    0   
[5,] 7.8      7.5  6.333333 5.5  6   

Por lo menos, podemos constatar que a medida que nuestra diagonal crece, la matriz de dimensión n-1 se sigue manteniendo.

No creo que sea más "performante" que tú solución, pero al menos es una aproximación distinta al problema.

1
  • Me agrada @Patricio Moracho el análisis que has planteado: "ver la matriz como un vector plano". Adicionalmente, si la diagonal tiene longitud 1 es soportado por vuestra función, algo que en la variante que planteé, simplemente no corre, a menos claro, que le agregue condicionales, pero eso era justo lo que se quería evitar. Commented el 15 mar. 2019 a las 20:13
1

En lugar de un 2 bucles explícitos anidados se pueden usar bucles implícitos con sapply(). En este caso creo que los for son más claros. Esta función es muy compacta y hace el trabajo, pero es bastante enredada.

promedio_matriz_triangular <- function(diagonal) {

  #Defino una función que crea cada columna de la matriz de salida
  #head() va reduciendo el largo de diagonal a medida que "bajo" una fila
  media_col <- function(x) sapply(seq_along(x), function(i) mean(head(x, i)))

  #Paso esa función sobre una versión cada vez más corta de diagonal
  #el rev() es necesario por algún motivo
  #el c() es para rellenar de 0 a la izquierda, igualando el largo de los vectores
  #sapply() coerciona a matriz por defecto
  sapply(rev(seq_along(diagonal)), function(j) c(rep(0, length(diagonal)-j), media_col(tail(diagonal, j))))
}

Eliminé la respuesta original por ser incorrecta (era para un problema diferente)

3
  • 1
    Hubert, El código está pensado para recibir directamente la matriz, no la diagonal. La función diag() regresa un resultado inesperado si le pasas un vector o matriz no rectangular. Más tarde lo corrijo con código defensivo.
    – mpaladino
    Commented el 15 mar. 2019 a las 16:12
  • Hubert, la respuesta anterior estaba equivocada, ya que no promediaba el intervalo entre los puntos de la diagonal. Esta es complicada y dudo que sea muy rápida, pero al menos da el resultado esperado.
    – mpaladino
    Commented el 16 mar. 2019 a las 1:44
  • Vale @mpaladino, lo interesante es mostrar diferentes aproximaciones al caso expuesto para que la comunidad se benefice y genial que uses tail en la solución, por lo general uno emplea esta función para ver si la data está completa, bien ahí. Commented el 16 mar. 2019 a las 3:30

Tu Respuesta

By clicking “Publica tu respuesta”, you agree to our terms of service and acknowledge you have read our privacy policy.

¿No es la respuesta que buscas? Examina otras preguntas con la etiqueta o formula tu propia pregunta.