[R] Trace of product of matrices

[Spencer Graves]

On 10/19/2014 8:42 AM, peter dalgaard wrote:
>> On 19 Oct 2014, at 16:43 , Wagner Bonat <wbonat at gmail.com> wrote:
>>
>> Dear,
>>
>> I have to compute the trace of a product between four matrices. For
>> example, I know the matrices Wi, Wj and C, I need to compute this
>>
>> -trace(Wi%*%C^-1%*%Wj%*%C^-1)
>>
>>
>> I would like to avoid compute the complete matrix and after take the
>> diagonal, something like
>>
>> sum(diag( solve(Wi,C)%*% solve(Wj,C)))
> <this can't be right: it is C that is the invertible matrix>
>
>> Any idea is welcome.
>>
> The usual "trick" is that the trace of a matrix product is the inner product in matrix space, which is just the sum of the  elementwise products
>
> tr(AB) = tr(BA) = sum_i sum_j a_ij b_ij.
>
> In R, this becomes simply sum(A*B) -- notice that the ordinary product is used, not %*%. So presumably, you are looking for
>
> sum(solve(C, Wi) * solve(C, Wj))

missing a transpose?


> A <- matrix(1:4, 2)
> sum(diag(A %*% A))
[1] 29
> sum(A*A)
[1] 30
> sum(diag(A %*% t(A)))
[1] 30
>

and  sum(A*t(A))
[1] 29


	  This is the answer you want, I think.  [Sorry:  I hit send before I 
thought through what I was doing.  Suppose A is a row vector and B a 
column vector.  Then the elementwise product A*B returns an error in R 
but A%*%B is a scalar while B%*%A is not, and sum(A*t(B)) = sum(t(A)*B) 
= trace(A%*%B).

When A and B are not square but the product is, then A*B = the 
elementwise product of A and B is not defined, but A*t(B) is, and so is 
t(A)*B.  When A


       See "Trace of a product" in the Wikipedia article on "Trace 
(linear algebra)" 
(https://en.wikipedia.org/wiki/Trace_%28linear_algebra%29#Trace_of_a_product). 



       Spencer

>> Thanks
>>
>> --
>> Wagner Hugo Bonat
>> LEG - Laboratório de Estatística e Geoinformação
>> UFPR - Universidade Federal do Paraná
>>
>> 	[[alternative HTML version deleted]]
>>
>> ______________________________________________
>> R-help at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.