如何在R中构建最优转移矩阵及用MDPtoolbox实现策略迭代

现有3组转移矩阵（P1、P2、P3）与奖励矩阵（R1、R2、R3），已计算得到各矩阵行乘积和的最大值：q1=8、q2=16、q3=7。以下是在R中构建目标转移矩阵P_NEW，以及使用MDPtoolbox完成策略迭代与改进的具体步骤。

初始化矩阵代码 #

# 策略1的初始变量设置
P1 <- matrix(c(1/2, 1/4, 1/4,
               1/16,3/4,3/16,
               1/4,1/8,5/8), nrow = 3, byrow = TRUE);rownames(P1)=c("1","2","3");P1
R1 <- matrix(c(10,4,8,
               8,2,4,
               4,6,4), nrow = 3, byrow = TRUE);rownames(R1)=c("1","2","3");R1


# 策略2的初始变量设置
P2 <- matrix(c(1/2, 0, 1/2,
               1/16,7/8,1/16), nrow = 2, byrow = TRUE);rownames(P2)=c("1","2");P2
R2 <- matrix(c(14,0,18,
               8,16,8), nrow = 2, byrow = TRUE);rownames(R2)=c("1","2");R2



# 策略3的初始变量设置
P3 <- matrix(c(1/4,1/4,1/2,
               1/8,3/4,1/8,
               3/4,1/16,3/18), nrow = 3, byrow = TRUE);rownames(P3)=c("1","2","3");P3
R3 <- matrix(c(10,2,8,
               6,4,2,
               4,0,8), nrow = 3, byrow = TRUE);rownames(R3)=c("1","2","3");R3

# 计算各矩阵行乘积和的最大值
q1 <- max(unname(rowSums(P1*R1)))
q2 <- max(unname(rowSums(P2*R2)))
q3 <- max(unname(rowSums(P3*R3)))

构建目标转移矩阵P_NEW #

1. 定位最大值对应的行 #

首先找到每个矩阵中，行乘积和等于对应最大值的行索引：

# 定位P1中rowSums(P1*R1)等于q1的行
row_p1 <- which(rowSums(P1*R1) == q1)
# 定位P2中rowSums(P2*R2)等于q2的行
row_p2 <- which(rowSums(P2*R2) == q2)
# 定位P3中rowSums(P3*R3)等于q3的行
row_p3 <- which(rowSums(P3*R3) == q3)

2. 提取行并构建P_NEW #

将找到的行按顺序组合，即可得到目标矩阵：

# 提取对应行并合并
selected_rows <- rbind(P1[row_p1, ], P2[row_p2, ], P3[row_p3, ])
# 转换为3行的矩阵
P_NEW <- matrix(selected_rows, nrow = 3, byrow = TRUE)

# 查看结果
P_NEW

运行后输出与目标矩阵一致：

[,1]  [,2]  [,3]
[1,] 0.50 0.25 0.25
[2,] 0.50 0.00 0.50
[3,] 0.25 0.25 0.50

使用MDPtoolbox进行策略迭代与改进 #

1. 安装并加载包 #

# 首次运行需安装包
install.packages("MDPtoolbox")
# 加载包
library(MDPtoolbox)

2. 定义MDP核心要素 #

设定状态数、动作数，构建转移概率数组和奖励数组：

# 状态数
S <- 3
# 动作数（此处先以P_NEW对应1个动作为例）
A <- 1

# 构建转移概率数组：维度为[状态数, 状态数, 动作数]
P <- array(0, dim = c(S, S, A))
P[, , 1] <- P_NEW

# 构建奖励数组：使用各状态的最大即时奖励（即q1、q2、q3）
R <- matrix(c(q1, q2, q3), nrow = S, ncol = A)

3. 执行策略迭代 #

调用mdp_policy_iteration函数完成策略迭代：

# 执行策略迭代（折扣因子gamma设为0.9，可根据需求调整）
pi_result <- mdp_policy_iteration(P, R, discount = 0.9)

# 查看输出结果
# 最优策略
cat("最优策略：", pi_result$policy, "\n")
# 最优值函数
cat("最优值函数：", pi_result$V, "\n")

4. 手动策略改进 #

若需手动执行策略改进，可基于当前值函数计算每个状态的最优动作：

gamma <- 0.9
# 计算每个状态-动作对的Q值
q_values <- matrix(0, nrow = S, ncol = A)
for (s in 1:S) {
  for (a in 1:A) {
    q_values[s, a] <- R[s, a] + gamma * sum(P[s, , a] * pi_result$V)
  }
}

# 选择每个状态下Q值最大的动作作为改进后的策略
improved_policy <- apply(q_values, 1, which.max)
cat("改进后的策略：", improved_policy, "\n")

内容的提问来源于stack exchange，提问作者Homer Jay Simpson