马尔可夫链用于程序生成

注意：该项目的完整源代码可以在[ 这里 ]找到。 由于它是一个较大项目的一部分，因此建议您在发布本文时或/source/helpers/arraymath.h文件以及/source/world/blueprint.cpp 。

在本文中，我想详细讨论在3D建筑物和其他系统的过程生成中使用马尔可夫链和统计的原理。

我将解释系统的数学基础，并尝试使解释尽可能笼统，以便您可以在其他情况下应用此概念，例如生成2D地牢。 解释中将附带图片和源代码。

该方法是用于满足特定要求的系统过程生成的一种通用方法，因此，我建议至少阅读第一节的末尾，以便您可以理解该技术是否对您的情况有用，因为下面将解释必要的要求。

结果将用在我的体素引擎中，以便任务机器人可以先建造建筑物，然后建造城市。 本文结尾处有一个示例！


带有结果的几个示例。

基础知识

马尔可夫链

马尔可夫链是系统沿其移动的状态序列，用时间的转变来描述。 状态之间的转换是随机的，也就是说，它们由概率来描述，而概率是系统的特征。

系统由状态空间定义，状态空间是所有可能的系统配置的空间。 如果正确地描述了系统，那么我们也可以将状态之间的转换描述为离散的步骤。

应该注意的是，从系统的一种状态来看，通常会有几种可能的离散过渡，每种过渡都会导致系统处于不同状态。

从状态i转换为状态j的概率等于：

$$

$$P_ {ij}$$

马尔可夫过程是借助转移给它的概率研究状态空间的过程。

重要的是，马尔可夫过程“没有记忆”。 这仅意味着从当前状态过渡到新状态的可能性仅取决于当前状态，而不再取决于先前访问的任何其他条件。

$$

$$P_ {ij} = P（i，j）$$

示例：文本生成

系统是位序列。 状态空间是所有可能的位序列。 离散转换将一位从0更改为1或从1更改为0。如果系统具有n位，则从每种状态我们都有n种可能的转换到新状态。 马尔可夫过程将包括通过使用某些概率更改序列中的位值来研究状态空间。

示例：天气预报

该系统是当前的天气状况。 状态空间是天气可能处于的所有可能状态（例如，“多雨”，“多云”，“晴天”等）。 过渡将是从任何状态到可以设置过渡可能性的另一种状态的转换（例如，“不管今天天气晴朗，明天是否下雨，无论昨天的天气如何？”的概率是多少？）。

马尔可夫链的蒙特卡罗方法

由于状态之间的转换由概率决定，因此我们还可以设置“稳定”处于任何状态的概率（或者，如果时间趋于无穷大，则我们处于特定状态的平均时间）。 这是状态的内部分布。

然后，用于马尔可夫链的蒙特卡洛算法（Markov-Chain Monte-Carlo，MCMC）是一种用于从状态空间获取样本的技术。 采样（sampling）是指在考虑内部分布的基础上，根据选择的概率来选择状态。

他们说，处于某种状态的可能性与某个成本函数成比例*，该成本函数给出了系统所处当前状态的“估计”。 可以认为，如果成本低，则处于该状态的可能性高，并且该比率是单调的。 成本函数定义如下：

$$

$$R（i）$$

注意：由于比例不一定是线性的，因此不确定是否正确使用“比例”一词。

然后，来自状态分布的样本将以较高的概率返回具有低成本（或良好的“估计”）的配置！

即使具有极大的状态空间（可能是无限的，但“可能是无限的”），不管系统的复杂程度如何，只要给我们足够的时间进行收敛，MCMC算法都将找到低成本的解决方案。

对状态空间的这种研究是随机优化的标准技术，并且在机器学习等领域有许多应用。

吉布斯分布

注意：如果您不清楚此部分，可以安全地跳过它。 您仍然可以利用系统的实施。

在确定了可能情况的成本之后，我们如何确定这种情况的可能性？

解决方案： Gibbs分布是给定约束集的最大熵分布。

从本质上讲，这意味着如果我们对系统的概率设置许多约束，则吉布斯分布将创建有关分布形状的“最少假设”。

注意：吉布斯分布也是对约束变化敏感度最低的分布（根据Kullback-Leibler散度度量标准）。

我们对状态分布施加的唯一限制是成本函数，因此我们在Gibbs分布中使用它来计算状态之间的转移概率：

$$

$$P_ {ij} = \ exp（-\ frac {R（j）-R（i）} {T}）* \ frac {1} {Z_i}$$

其中Z是从状态i开始的过渡集合的分区函数。 这是一个标准化因子，可确保任何状态的转移概率之和为1。

$$

$$Z_i = \ sum_j（P_ {ij}）$$

请注意，如果我们确定下一个状态将是相同状态，则相对成本为零，即归一化后的概率将为非零（由于具有指标的分布形状）！ 这意味着在许多转换中，必须包括状态不变的概率。

还值得注意的是，吉布斯分布由“计算温度” T参数化。

在状态空间研究中使用概率的关键优势之一是，系统可以执行到更昂贵状态的转换（因为它们具有非零转换概率），从而将算法转变为“非贪婪”优化方法。

请注意，随着温度趋于无穷大，任何单个跃迁的可能性都趋向于以这样一种方式统一：当状态的所有跃迁的概率集被归一化时，尽管它们的成本更高，但它们变得同等概率（或吉布斯分布接近正态分布）！

随着计算温度趋近于零，具有较低成本的转换变得更有可能，即，首选转换的可能性增加。

在进行状态空间的研究/优化时，我们逐渐降低温度。 该过程称为“模拟退火” 。 因此，一开始我们可以轻松摆脱本地最小值，最后我们可以选择最佳解决方案。

当温度足够低时，除无过渡的可能性外，所有概率都趋于零！

这是因为仅不存在过渡具有零成本差异，即处于相同状态不依赖于温度。 由于在T = 0处的指数函数的形状，这成为具有非零值的唯一概率，也就是说，在归一化之后，它变成了单位。 因此，我们的系统将收敛到稳定点，不再需要进一步的冷却。 这是使用吉布斯分布的概率生成的不可或缺的特性。

可以通过更改冷却速度来调整系统的收敛过程！

如果冷却速度较慢，那么结果是我们通常会以较低的成本（在一定程度上）得出解决方案，但要付出更多的收敛步骤。 如果冷却速度更快，那么早期的系统更有可能陷入成本较高的次区域的陷阱，也就是说，我们将获得“次优”的结果。

因此，马尔可夫过程不仅会产生随机结果，而且还会尝试产生“良好”结果，并且很有可能成功！

根据任意成本函数的定义，不必存在唯一的最优值。 这种概率优化方法仅生成接近最佳值的方法，试图使成本函数最小化，并且由于抽样，每次生成的结果都不同（如果随机数生成器的种子不同）。

可以使用逆变换方法对离散过渡集的质量分布函数执行采样过程本身。 稍后，我将展示如何完成此操作。

程序生成

此方法对程序生成有何用处？

在某些系统中，通常很难定义产生良好结果的简单算法，尤其是在复杂系统的情况下。

设置任意的生成规则不仅困难，而且仅受我们的想象力和边界案例处理的限制。

如果系统满足一组特定的要求，那么MCMC的使用将使我们不必担心算法或规则的选择。 相反，我们定义了一种生成任何可能结果的方法，并根据其“评估”自觉地选择一个好的结果。

提出以下要求：

• 系统可以处于离散（可能是无限）状态配置。
• 我们可以定义状态之间的离散过渡。
• 我们可以设置一个成本函数来估计系统的当前状态。

在我下面，我将给出一些其他示例，在我看来，可以应用此方法。

实作

伪代码

在我们的实现中，我们希望实现以下目标：

• 设置系统状态。
• 将所有可能的过渡设置为下一个状态。
• 计算当前状态的成本。
• 计算所有可能的下一个状态（或其子集）的成本。
• 使用吉布斯分布，计算转变的概率。
• 使用概率进行样本（样本）转换。
• 执行采样的过渡。
• 降低计算温度。
• 重复步骤，直到获得满意的结果。

MCMC算法以伪代码的形式如下：

 // MCMC    T = 200; //     State s = initialState(); Transitions t[n] = {...} //n   thresh = 0.01; //  ( ) // ,      while(T > thresh){ //    curcost = costfunc(s); newcost[n] = {0}; // newcost   0 probability[n] = {0}; //     0 //  for(i = 0; i < n; i++){ newcost[i] = costfunc(doTransition(s, t[i])); probability[i] = exp(-(newcost[i] - curcost)/T); } //  probability /= sum(probability); //  sampled = sample_transition(t, probability); //  s = doTransition(s, sampled); //  T *= 0.975; } 

3D建筑物生成

状态空间与过渡

为了生成3D建筑物，我生成了许多具有边界框指定体积的房间。

 struct Volume{ //   glm::vec3 a; glm::vec3 b; void translate(glm::vec3 shift); int getVol(); }; //   int Volume::getVol(){ return abs((bx-ax)*(by-ay)*(bz-az)); } //    void Volume::translate(glm::vec3 shift){ a += shift; b += shift; } 

如果我生成n个房间，则系统的状态是3D空间中边框的配置。

应该注意的是，这些卷的可能配置是无止境的，但无穷无尽的（它们可以在无限长的时间内列出）！

 //  (  !) std::vector<Volume> rooms; // N  for(int i = 0; i < n; i++){ //  Volume x; xa = glm::vec3(0); xb = glm::vec3(rand()%4+5); //   //  . rooms.push_back(x); } //... 

许多可能的转换是将房间在六个空间方向之一上移动一步，包括不进行转换：

 //... //   std::array<glm::vec3, 7> moves = { glm::vec3( 0, 0, 0), //   ! glm::vec3( 1, 0, 0), glm::vec3(-1, 0, 0), glm::vec3( 0, 1, 0), glm::vec3( 0,-1, 0), glm::vec3( 0, 0, 1), glm::vec3( 0, 0,-1), }; //... 

注意：重要的是要使系统能够保持其当前状态！

成本函数

我希望3D空间中的体积表现得像“磁铁”，即：

• 如果它们的体积相交，那么这是不好的。
• 如果它们的表面接触，则很好。
• 如果他们根本不触摸，那就不好了。
• 如果他们碰到地板，那很好。

对于3D空间中的两个长方体，我们可以轻松地定义一个边界框：

 Volume boundingBox(Volume v1, Volume v2){ //   Volume bb; //   bb.ax = (v1.ax < v2.ax)?v1.ax:v2.ax; bb.ay = (v1.ay < v2.ay)?v1.ay:v2.ay; bb.az = (v1.az < v2.az)?v1.az:v2.az; //   bb.bx = (v1.bx > v2.bx)?v1.bx:v2.bx; bb.by = (v1.by > v2.by)?v1.by:v2.by; bb.bz = (v1.bz > v2.bz)?v1.bz:v2.bz; return bb; } 

使用体积的边界框，我们可以计算一个3D向量，该向量为我提供有关两个体积的交点的信息。

如果平行六面体沿一侧的长度大于沿该侧的两个体积的长度之和，则它们从该侧不接触。 如果它们相等，则曲面相接触；如果它们相等，则曲面相交。

 //    3  glm::vec3 overlapVolumes(Volume v1, Volume v2){ //      Volume bb = boundingBox(v1, v2); //  glm::vec3 ext1 = glm::abs(v1.b - v1.a); // v1  3  glm::vec3 ext2 = glm::abs(v2.b - v2.a); // v2  3  glm::vec3 extbb = glm::abs(bb.b - bb.a); //   //  return ext1 + ext2 - extbb; } 

该代码用于计算我形成加权量的数量，最终将其用作成本。

 int volumeCostFunction(std::vector<Volume> rooms){ // int metric[6] = { 0, //   0, //   0, //     0, //     0, // ,   0};//    int weight[6] = {2, 4, -5, -5, -5, 5}; //    for(unsigned int i = 0; i < rooms.size(); i++){ //     for(unsigned int j = 0; j < rooms.size(); j++){ //    ,  . if(i == j) continue; //    . glm::vec3 overlap = overlapVolumes(rooms[i], rooms[j]); //   glm::vec3 posOverlap = glm::clamp(overlap, glm::vec3(0), overlap); metric[0] += glm::abs(posOverlap.x*posOverlap.y*posOverlap.z); //   //   glm::vec3 negOverlap = glm::clamp(overlap, overlap, glm::vec3(0)); metric[1] += glm::abs(negOverlap.x*negOverlap.y*negOverlap.z); //   //   if(overlap.y == 0){ metric[2] += overlap.x*overlap.z; } //   (X) if(overlap.x == 0){ //      0,   , .. overlap.z = 0 metric[3] += overlap.z*overlap.y; } //   (Z) if(overlap.z == 0){ //      0,   , .. overlap.x = 0 metric[4] += overlap.x*overlap.y; } } //  ,   -   . if(rooms[i].ay == 0){ //  ,  ,    . metric[4] += rooms[i].ax*rooms[i].az; } //,     ! if(rooms[i].ay < 0){ //,        if(rooms[i].by < 0){ metric[5] += rooms[i].getVol(); } else{ metric[5] += abs(rooms[i].ay)*(rooms[i].bz-rooms[i].az)*(rooms[i].bx-rooms[i].ax); } } } // Metric * Weights return metric[0]*weight[0] + metric[1]*weight[1] + metric[2]*weight[2] + metric[3]*weight[3] + metric[4]*weight[4] + metric[5]*weight[5]; } 

注意：“正相交体积”是指实际相交的体积。 “负相交体积”表示它们根本不接触，并且相交由位于3D空间中两个长方体的两个最近点之间的空间中的体积定义。

权重的选择应优先考虑一件事，而对另一件事要罚款。 例如，在这里，我对位于地板下的房间进行了严格的罚款，同时还提高了其表面接触的那些建筑物的优先级（比对体积的交点进行罚款还要多）。

我为所有成对的房间生成成本，而忽略与自己配对的房间。

寻找低成本解决方案

如伪代码中所述执行收敛。 进行过渡时，我一次只能移动一个房间。 这意味着对于n个房间和7个可能的过渡，我需要计算7 * n个概率，然后从所有这些中选择，仅移动那个房间，这可能是最可取的。

  //  float T = 250; while(T > 0.1){ //   std::vector<std::array<double, moves.size()>> probabilities; //   () int curEnergy = volumeCostFunction(rooms); //      ... for(int i = 0; i < n; i++){ //    std::array<double, moves.size()> probability; //      ,     for(unsigned int m = 0; m < moves.size(); m++){ //        . rooms[i].translate(moves[m]); //      ! probability[m] = exp(-(double)(volumeCostFunction(rooms) - curEnergy)/T); //   rooms[i].translate(-moves[m]); } //       probabilities.push_back(probability); } //  ( ) double Z = 0; for(unsigned int i = 0; i < probabilities.size(); i++){ for(unsigned int j = 0; j < probabilities[i].size(); j++){ Z += probabilities[i][j]; } } //  for(unsigned int i = 0; i < probabilities.size(); i++){ for(unsigned int j = 0; j < probabilities[i].size(); j++){ probabilities[i][j] /= Z; } } //    (CDF) ( ) std::vector<double> cdf; for(unsigned int i = 0; i < probabilities.size(); i++){ for(unsigned int j = 0; j < probabilities[i].size(); j++){ if(cdf.empty()) cdf.push_back(probabilities[i][j]); else cdf.push_back(probabilities[i][j] + cdf.back()); } } //      std::random_device rd; std::mt19937 e2(rd()); std::uniform_real_distribution<> dist(0, 1); double uniform = dist(e2); int sampled_index = 0; //   CDF for(unsigned int i = 0; i < cdf.size(); i++){ //    ,   ... if(cdf[i] > uniform){ sampled_index = i; break; } } //     int _roomindex = sampled_index/moves.size(); int _moveindex = sampled_index%moves.size(); //  rooms[_roomindex].translate(moves[_moveindex]); // T T *= 0.99; // !!! } //!! //... 

; « ».

, (cumulative distribution function, CDF). 0 1. CDF , , « CDF », . :


代替连续功能，可能会有不连续的步骤。更多细节可以在这里阅读。

另外，我在3D空间中有房间体积数据！

我使用它们使用blueprint类将主题应用于已知的批量数据，以生成“模式”。这样房屋就可以出现了。先前的文章[ here ]（Habré的翻译）中描述了蓝图类。有关从这些卷中完整生成房屋的信息，请参见源代码。

结果

这种广义方法的结果非常好。我唯一需要设置的是成本函数中正确的优先级和惩罚权重。


使用此算法的建筑生成的几个示例以及适用于它们的主题。


( ).

, (3-5), .

, 3D- , MCMC.

, , , . , , .

. ( ).

, 3D- 2D-, .

, 2D- — , .

, , , , , .

MCMC, , ( , ..).

:

• ; , / .
• , , , .
• , !

: MCMC , « », NP- . — , , — !

Task-

, task- .

, , :


( , ). . ( ) . , . , . - , , .