R代码优化：几何布朗运动三叉树节点扩展至20层的实现方案

问题描述 #

我需要生成一棵展示标的资产遵循几何布朗运动演化的三叉树，目前的核心问题是节点数量以3^{n的规模增长，当n=4时就得写27个循环，效率极低。我希望优化现有的R代码，能把节点层数扩展到20层（对应矩阵行数为3}20），同时保留当前的树结构，不能改成常见的涨跌平结构。以下是我的现有代码：

sigma=0.15  # questo l'ho scelto io 
hist_sigma<-sd(AAPL$AAPL.Close)
S<-matrix(0,nrow = 27,ncol = 4)
S0=S[1,1]<-last(AAPL$AAPL.Close)
T=90/360
set.seed(22)
Z<-matrix(data = 0,nrow = 27, ncol = 3)
for (i in 1:27) {
  Z[i,]<-rnorm(3,mean = 0,sd=1)
  
}



for (i in 1:3) {
  for (j in 2:2) {
    S[i,j]<-S[1,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
  }
}

 for (i in 1:3) {
  for (j in 3:3) {
    S[i,j]<-S[1,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
    
  }
 }

for (i in 4:6) {
  for (j in 3:3) {
    S[i,j]<-S[2,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
    
  }
}

for (i in 7:9) {
  for (j in 3:3) {
    S[i,j]<-S[3,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
    
  }
}

for (i in 1:3) {
  for (j in 4:4) {
    S[i,j]<-S[1,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
    
  }
}
for (i in 4:6) {
  for (j in 4:4) {
    S[i,j]<-S[2,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
    
  }
}

for (i in 7:9) {
  for (j in 4:4) {
    S[i,j]<-S[3,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
    
  }
}
for (i in 10:12) {
  for (j in 4:4) {
    S[i,j]<-S[4,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
    
  }
}
for (i in 13:15) {
  for (j in 4:4) {
    S[i,j]<-S[5,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
    
  }
}

for (i in 16:18) {
  for (j in 4:4) {
    S[i,j]<-S[6,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
    
  }
}
for (i in 19:21) {
  for (j in 4:4) {
    S[i,j]<-S[7,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
    
  }
}

for (i in 22:24) {
  for (j in 4:4) {
    S[i,j]<-S[8,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
    
  }
}

for (i in 25:27) {
  for (j in 4:4) {
    S[i,j]<-S[9,j-1]*(exp((r-1/2*sigma^2)*(120/360) + sigma*(sqrt((120/360))*Z[i,j-1])))
    
  }
}

优化思路 #

1. 利用三叉树层级索引规律：每一层节点数是上一层的3倍，第j-1层的第i个节点对应第j层的3*(i-1)+1到3*i行，通过这个数学关系自动映射父节点与子节点，无需手动编写循环。
2. 向量化运算提升效率：R的向量运算比嵌套循环快得多，批量处理每个父节点对应的3个子节点。
3. 预分配内存：一次性计算总节点数并分配矩阵内存，避免动态扩容的性能损耗。

优化后的代码 #

# 参数配置
sigma <- 0.15
hist_sigma <- sd(AAPL$AAPL.Close)
S0 <- last(AAPL$AAPL.Close)
total_time <- 90/360  # 总时间跨度
n_layers <- 20        # 目标层数
dt <- total_time / n_layers  # 每层时间步长
r <- 0.05             # 补充原代码缺失的无风险利率，可按需调整

# 计算总节点数：3^(n_layers-1)，对应n_layers层的三叉树
total_nodes <- 3^(n_layers - 1)
# 预分配价格矩阵
S <- matrix(0, nrow = total_nodes, ncol = n_layers)
S[1, 1] <- S0  # 初始节点价格

# 生成所有随机数：每个节点需要(n_layers-1)个随机值（从第2层到第n层）
set.seed(22)
Z <- matrix(rnorm(total_nodes * (n_layers - 1), mean = 0, sd = 1),
            nrow = total_nodes, ncol = n_layers - 1)

# 预计算指数项的固定部分，减少重复计算
exp_const <- (r - 0.5 * sigma^2) * dt
exp_sigma <- sigma * sqrt(dt)

# 逐层生成节点价格
for (j in 2:n_layers) {
  prev_node_count <- 3^(j - 2)  # 上一层的节点数量
  for (i in 1:prev_node_count) {
    parent_price <- S[i, j-1]
    # 计算当前父节点对应的子节点行范围
    child_row_range <- (3*(i-1)+1):(3*i)
    # 批量计算子节点价格
    S[child_row_range, j] <- parent_price * exp(exp_const + exp_sigma * Z[child_row_range, j-1])
  }
}

关键说明 #

• 索引映射：通过3*(i-1)+1到3*i的行范围，精准对应每个父节点的3个子节点，完全保留原树结构。
• 性能优化：向量化计算子节点价格，避免了原代码中大量重复的循环逻辑，即使到20层也只需2层嵌套循环。
• 内存注意：3^19=1162261467个节点，每个数值占约8字节，总内存需求约9GB，需确保机器有足够内存。若内存不足，可考虑分块处理或使用稀疏矩阵（但会增加复杂度）。
• 时间步长修正：原代码中总时间90/360但每层时间步长120/360存在逻辑矛盾，优化后将总时间平均分配到每层，若需保持原时间步长，可将dt改为120/360，同时调整层数与总时间的匹配关系。

内容的提问来源于stack exchange，提问作者James