Llopis Llopis - 3 months ago 11
R Question

Sum of products which indice is not i an j

I am implementing the following function in the R code:

formula to be converted into code

So far I used:

sig.TOM <- function(adj, sig.adj) {
out <- matrix(nrow = nrow(adj), ncol = ncol(adj))
for (i in 1:nrow(adj)) {
for (j in 1:ncol(adj)) {
out[i,j] <- abs(adj[i, j] + sum(sig.adj[i, -c(i, j)]*sig.adj[-c(i, j), i]))/(
min(sum(sig.adj[-i, i]), sum(sig.adj[-j, j])) + 1 - abs(adj[i,j]))
}
}
return(out)
}


where ~a is the following mock matrix:

sig.adj <- structure(c(1, -0.418913311940584, 1, 0.947013383275973, -1,
-0.418913311940584, 1, -0.207962861914701, 0.584386281408348,
-0.687223049826016, 1, -0.207962861914701, 1, 0.763551721347657,
-0.0327147711077901, 0.947013383275973, 0.584386281408348, 0.763551721347657,
1, 0.284466543760789, -1, -0.687223049826016, -0.0327147711077901,
0.284466543760789, 1), .Dim = c(5L, 5L))


and
adj <- abs(sig.adj)
, where adj in the formula is described as a and sig.adj as ~a.

But the result is not symmetric as expected so I must have implemented it wrong, I have doubts on the summation part.

How can that sum of products of values when the indices are not i or j be implemented?

Answer

If I understood the formula you want:

sig.adj <- structure(c(1, -0.418913311940584, 1, 0.947013383275973, -1, 
                       -0.418913311940584, 1, -0.207962861914701, 0.584386281408348, 
                       -0.687223049826016, 1, -0.207962861914701, 1, 0.763551721347657, 
                       -0.0327147711077901, 0.947013383275973, 0.584386281408348, 0.763551721347657, 
                       1, 0.284466543760789, -1, -0.687223049826016, -0.0327147711077901, 
                       0.284466543760789, 1), .Dim = c(5L, 5L))
spec.mult <- function(A,B) {
  rA <- nrow(A); cB <- ncol(B)
  C <- A %*% B
  for (i in 1:rA) for (j in 1:cB) 
    C[i,j] <- C[i,j] - A[i,i]*B[i,j] - A[i,j]*B[j,j] + ifelse(i==j, A[i,i]*B[j,j], 0) 
  C
}

sig.TOM <- function(sig.adj) {
  adj <- abs(sig.adj)
  k <- colSums(adj) - abs(diag(adj))
  abs(adj + spec.mult(sig.adj, sig.adj)) / (outer(k,k, pmin) +1 - abs(adj))
}
sig.TOM(sig.adj)

p.s.: I copied the part ... +1 - abs(adj) from your question because I don't know if you want ... +1 - adj or ... +1 - sig.adj

Here is another implementation for symmetric matrix sig.adj (as in your case):

spec.mult <- function(A) {
  dA.A <- diag(A)*A
  crossprod(A) - dA.A - t(dA.A) + diag(diag(A)^2)
}

sig.TOM <- function(sig.adj) {
  adj <- abs(sig.adj)
  k <- colSums(adj) - abs(diag(adj))
  abs(adj + spec.mult(sig.adj)) / (outer(k,k, pmin) +1 - abs(adj))
}
sig.TOM(sig.adj)