0

Good day, Please i need assistance trying to generalize this code for any matrix the data PoliticalDemocracy is available in the lavaan library

What I tried

library(lavaan)
R<-chol(cor(PoliticalDemocracy[1:8]))
diag(R)<-NA
R<-round(R,3)
R
R[R==0]<-NA
R

d<-matrix(NA,nrow(R),ncol(R))
d[,1]<-acos(R[1,])
d


d[,2]<-acos(R[2,]/sin(d[,1]))
d[,3]<-acos(R[3,]/sin(d[,1])*sin(d[,2]))
d[,4]<-acos(R[4,]/sin(d[,1])*sin(d[,2])*sin(d[,3]))
d[,5]<-acos(R[5,]/sin(d[,1])*sin(d[,2])*sin(d[,3])*sin(d[,4]))
d[,6]<-acos(R[6,]/sin(d[,1])*sin(d[,2])*sin(d[,3])*sin(d[,4])*sin(d[,5]))
d[,7]<-acos(R[7,]/sin(d[,1])*sin(d[,2])*sin(d[,3])*sin(d[,4])*sin(d[,5])*sin(d[,6]))
d[,8]<-acos(R[8,]/sin(d[,1])*sin(d[,2])*sin(d[,3])*sin(d[,4])*sin(d[,5])*sin(d[,6])*sin(d[,7]))

I tried using the loop but was having errors and an empty matrix d


d<-matrix(NA,nrow(R),ncol(R))
d[,1]<-acos(R[1,])
d


for (i in 2:nrow(d)){
  d[,i]<-acos(R[i,])/prod(sin(d[,i-1]))
  d
}
d

@AllanCameron

  • Is its `acos(a/sin(b)*sin(c))` or is it `acos(a/(sin(b)*sin(c)))`? Check the paranthesis. – Onyambu Apr 20 '22 at 14:07
  • Also note that you have `NA` in every row. Meaning whenever you multiply, you will get NA. so what exactly are you trying to achieve?? – Onyambu Apr 20 '22 at 14:22
  • The lines above where I feel the matrix d, column by column is what I need, I used the NA since it's a triangular matrix and I wouldn't want the arccosine of 0 to return any value as it is 1.57, so trying to reduce the codes to a generalization so I can use it for any matrix irrespective of the length is the challenge now. – Etini Akpayang Apr 21 '22 at 06:13
  • You are multiplying by NA. Also you never responded as to whether it is `acos(a/sin(b)*sin(c))` or `acos(a/(sin(b)*sin(c)))` – Onyambu Apr 21 '22 at 06:15
  • Should be the second one, everything in a column in R, divided by the product of the sines, and then the arc cosine of the results. – Etini Akpayang Apr 21 '22 at 06:16
  • Then your whole code above is incorrect since there is no `paranthesis to enclose the divisor – Onyambu Apr 21 '22 at 06:18
  • Ow! That's right! I'll try that and get back to you much later as am not currently near a PC. – Etini Akpayang Apr 21 '22 at 06:19

1 Answers1

0

You could use Reduce function as shown below

(Reduce(function(x, y) 
  cbind(x, acos(y/matrixStats::rowProds(sin(x)))),
  data.frame(t(R)), init = matrix(rep(pi/2), nrow(R)))[,-1])
Onyambu
  • 67,392
  • 3
  • 24
  • 53