問題描述
如何在循環和 if 語句中使用遞歸公式 (How to use recursive formula in loop and if statement)
示例矩陣的偽代碼</em></p>
a = lower triangular matrix/dataframe 20x20
b = a 1x20 matrix/vector
c = the previous row result of the formula (recursive bit)
上我目前的代碼,
先創建變量如下圖:
Matrix C
1 | 2 |3 |4 |5 |6 |
1|i=j=1 |i<j=0 |...|...|...|...|
2|i>j=a(2,1)b(1)c(1,1)/b(2)‑b(1) |i=j=1 |...|...|...|...|
3|i>j= (a(3,1)b(1)c(1,1) + a(3,2)b(2)c(2,1))/(b(3)‑b(1) |i>j= a(3,2)b(2)c(2,2)/(b(3)‑b(2) |...|...|...|...|
4|i>j = a(4,1)b(1)c(1,1)+a(4,2)b(2)c(2,1)+a(4,3)b(3)c(3,1))/b(4)‑b(1) |i>j= a(4,2)b(2)c(2,2) + a(4,3)b(3)c(3,2)/b(4)‑b(2) |...|...|...|...|
5|... |... |...|...|...|...|
6|... |... |...|...|...|...|
</code></pre>
計算公式結果,我認為由於 R 的矢量化而有效
my_int <‑ 20
nr <‑ as.integer(my_int)
#create a n x n matrix with zeroes
a <‑ matrix(0, nr, nr)
# For each row and for each column, assign values based on position
# These values are the product of two indexes
for(i in 1:dim(a)[1]) {
for(j in 1:dim(a)[2]) {
a[i,j] = if(i<j) {
0
}else if(i==j) {
1
}else {
3
}
}
}
# make into dataframe
mymat <‑ data.frame(mymat)
# create b variable
z <‑ rep(1:20)
# made it into lower diagonal matrix to make it easier to work with
b <‑ matrix(0, length(z), length(z))
b[lower.tri(b, diag = TRUE)] <‑ z[sequence(length(z):1)]
b
# create matrix for formula to operate with
# create variable "c"
c <‑ matrix(0, nr, nr)
# For each row and for each column, assign values based on position
# These values are the product of two indexes
for(i in 1:dim(c)[1]) {
for(j in 1:dim(c)[2]) {
mymat2[i,j] = if(i<j) {
0
}else if(i==j) {
1
}else {
5 # place holder for now
}
}
}
</code></pre>
然後我的問題是如何將其插入到 if 語句中以創建遞歸定義的矩陣,如下所示
sum(ablag(c), na.rm = TRUE)/(b[,j]‑b[I,])
</code></pre>
這給出錯誤
mymat2[i, j] <‑ if (i
如果有幫助,我可以鏈接一個 Excel 電子表格,該公式適用於該電子表格。它只能通過大量手動輸入等在 excel 中實現。
使用真實數據的預期結果示例
# calculate recursively defined lower triangular matrix
c <‑ matrix(0, nr, nr)
# For each row and for each column, assign values based on position
# These values are the product of two indexes
for(i in 1:dim(c)[1]) {
for(j in 1:dim(c)[2]) {
mymat2[i,j] = if(i<j) {
0
}else if(i==j) {
1
}else {
sum(a*b*lag(c), na.rm = TRUE)/(b[,j]‑b[I,]) # formula for calculation of values for lower triangular matrix
}
}
}
</code></pre>
我如何計算“c”的示例 在 excel 中
a = 5x5 lower triangle matrix
0|0 |0 |0 |0
5|0 |0 |0 |0
5|0.56|0 |0 |0
5|0.20|0.61|0 |0
5|0.06|0.16|0.61|0
b = 1x5 matrix/vector
0.27917|0.499|0.83|1.191|1.48
c = recursive matrix results
1 |0 |0 |0 |0
6.36|1 |0 |0 |0
5.77|0.84|1 |0 |0
5.50|0.77|1.43|1 |0
5.3 |0.72|1.80|2.46|1
</code></pre>
上面的代碼顯示瞭如何在 excel 中通過索引和手動設置大量單元格位置來計算它。
參考解法
方法 1:
Here's a possibly dynamic programming(?) solution. I've only verified a few 'cells' so you'll need to verify that it's giving results you expected.
# FUNCTIONS
vpad <‑ function(vec) {
extend <‑ 20 ‑ length(vec)
c(vec, rep(0, extend))
}
clip <‑ function(x) { x[x >= 0 & x <= 20] }
pad_a <‑ function(a, iter) {
mat <‑ cbind(a[, iter:20], matrix(0, 20, iter‑1))
mat * lower.tri(mat)
}
pad_b <‑ function(b, iter) {
tmp <‑ vpad(b[iter:20])
mat <‑ matrix(rep(tmp, each=20), 20, 20)
mat * lower.tri(mat)
}
pad_c <‑ function(c, iter) {
L <‑ 21 ‑ iter
i <‑ clip(seq(L) + iter ‑ 1)
j <‑ seq(i)
v <‑ vpad(mapply(function(i, j) c[i, j], i, j))
mat <‑ matrix(rep(v, each=20), 20, 20)
mat <‑ mat * (lower.tri(mat) + diag(20))
mat
end <‑ 21 ‑ iter
mat1 <‑ rbind(matrix(0, iter ‑ 1, 20), mat[1:end, ])
mat1 * lower.tri(mat1)
}
# MAKE FAKE DATA
set.seed(1)
m <‑ matrix(sample(400)/400, 20, 20)
a <‑ m * lower.tri(m)
b <‑ sample(20)/21
# INITIALIZE VARS
add_prev <‑ matrix(rep(0, 200), 20, 20)
init_mat <‑ lower.tri(m) + diag(20)
overwrite_mat <‑ init_mat
for (i in 1:20) {
overwrite_mat <‑ (pad_a(a, i) * pad_b(b, i) * pad_c(overwrite_mat, i)) + add_prev
add_prev <‑ overwrite_mat
}
denom1 <‑ matrix(rep(b, times=20), 20, 20)
denom1 <‑ denom1 * lower.tri(denom1) + diag(20)
denom2 <‑ matrix(rep(b, each=20), 20, 20)
denom2 <‑ denom2 * lower.tri(denom2) + diag(20)
final <‑ ((overwrite_mat / denom1) ‑ denom2) + diag(20)
final <‑ replace(final, is.nan(final), 0)
Output (just showing upper left (5,5) matrix)
[,1] [,2] [,3] [,4] [,5]
[1,] 0.0000000 0.00000000 0.00000000 0.00000000
[2,] ‑0.5129073 0.00000000 0.00000000 0.00000000
[3,] 0.6738060 ‑0.05926190 0.00000000 0.00000000
[4,] 9.2856557 5.97660268 ‑0.15059524 0.00000000
[5,] 0.1971807 0.01591345 ‑0.14775033 ‑0.08718254
(by Biggerthanpenny、CPak)
參考文件