2017-10-20 32 views
0

j'ai une norme qui est décrite par la matrice sigmamultiplication norme/matrice vectoriced

sigma <- matrix(c(1,0.5,0,0.5,1,0,0,0,1),3,3)) 

Pour calculer la norme d'un vecteur I COMPUTE

t(x) %*% sigma %*% x 

qui fonctionne très bien pour un vecteur, par exemple x = 1:3.

Cependant, je veux calculer la norme de nombreux vecteurs en même temps, ce qui est je

x <- t(matrix(rep(1:3, 10),3,10)) 

(bien sûr rempli d'entrées différentes).

Existe-t-il un moyen de calculer la norme de chaque vecteur simultanément? I.e. quelque chose comme

lapply(1:10, function(i) t(x[i,]) %*% sigma %*% x[i,]) 

Répondre

3

Vous pouvez faire:

sigma <- matrix(c(1,0.5,0,0.5,1,0,0,0,1),3,3) 
x <- t(matrix(rep(1:3, 10),3,10)) 

mynorm <- function(x, sig) t(x) %*% sig %*% x 
apply(x, 1, mynorm, sig=sigma) 

Voici une variante avec tcrossprod():

mynorm <- function(x, sig) tcrossprod(x, sig) %*% x 
apply(x, 1, mynorm, sig=sigma) 

Et est la référence ici (y compris variantes de la solution de compute only diagonals of matrix multiplication in R, merci à @Benjamin pour le lien):

mynorm1 <- function(x, sig) t(x) %*% sig %*% x 
mynorm2 <- function(x, sig) tcrossprod(x, sig) %*% x 

microbenchmark(n1=apply(x, 1, mynorm1, sig=sigma), 
       n2=apply(x, 1, mynorm2, sig=sigma), 
       n3 = colSums(t(x) * (sigma %*% t(x))), 
       n4 = rowSums(x * t(sigma %*% t(x))), 
       n5 = rowSums(x * (x %*% t(sigma))), 
       n6 = rowSums(x * tcrossprod(x, sigma)), 
       Eugen1 = diag(x %*% sigma %*% t(x)), 
       Eugen2 = diag(x %*% tcrossprod(sigma, x)), 
       unit="relative") 
2

Que pensez-vous de cette simple multiplication matricielle:

diag(t(x) %*% sigma %*% x) 

Edit: après la matrice multiplications vous avez besoin de la diagonale (bien sûr).

Ensuite, il est plus rapide que les solutions avec appliquent

+0

Cela me donne "Erreur dans t (x)% *% sigma: arguments non concordantes" –

+0

Ok, alors x% *% sigma% *% t (X). Cela a fonctionné pour moi. – EugenR

+0

@EugenR Vous calculez des éléments inutiles. – jogo

2

Cela devrait faire

> sigma <- matrix(c(1,0.5,0,0.5,1,0,0,0,1),3,3) 
> x <- t(matrix(rep(1:30, 10),3,10)) 
> 
> # should give 
> t(x[1, ]) %*% sigma %*% x[1, ] 
    [,1] 
[1,] 16 
> t(x[2, ]) %*% sigma %*% x[2, ] 
    [,1] 
[1,] 97 
> 
> # which you can get by 
> rowSums((x %*% sigma) * x) 
[1] 16 97 250 475 772 1141 1582 2095 2680 3337 
+1

La question est proche d'être un doublon de https://stackoverflow.com/q/21708489/5861244 –