重要性采样(Importance Sampling)简介和简单样例实现

在渲染领域，重要性采样这个术语是很常见的，但它究竟是什么呢？我们首先考虑这样的一种情况：

如果场景里有一点P，我们想计算P点的最终颜色，根据全局照明的概念，P点的颜色是由所有投射到P点的所有光线所影响的。但是我们在利用蒙特卡洛积分去计算投射光线的总和会出现一个问题，如上图场景里有两盏灯，当蒙特卡洛采样较少的情况，可能导致灯光方向未被采样，这样计算的结果会和实际产生巨大的偏差，尤其是在采样函数(f(X))有突峰的情况。如下图所示：

而在采样函数是常数时，则完全没有该问题：

当然，我们几乎不会使用是常数的采样函数。可是我们可以这样计算蒙特卡洛积分：

我们寻找一个pdf(X)函数是f(X)的倍数，这样采样函数会变成如下图所示的常数了：

该图正好是我们上面两盏灯场景的光线强度采样函数，当我们取 ${f}'(X)$ 为 ${f}(X)$ 的倍数，则我们的采样就能像采样常数那样进行。当我们的pdf选择上图中所示的函数时，我们的X则正好有在两盏灯光附近采样较多样本，其他地方采样较少的样本的特点。这就是我们所说的重要性采样。虽然这不是均匀采样，但是我们在 ${f}(X)$ 大的时候，我们的pdf值也大， ${f}(X)$ 小的时候，pdf值也小，类似相互补足了，也导致我们的函数最终能够正确的收敛，而且具有更小的方差。但实际上很难找到完全是倍数的pdf，但我们可以找形状十分类似的pdf。

我们现在用一个例子：

我们使用蒙特卡洛积分和重要性采样来计算该积分。虽然该积分可以很快速的用计算获得：

但是我们的目的是为了验证重要性采样的快速收敛性。我们选择两个pdf进行对比，如下图所示：

一个是 $p(x)=1/\pi$ 的均匀分布，另一个则是 ${p}'(x)=8x/\pi^{2}$ 形状与 $f(x)$ 类似的分布。在我们用代码验证之前，我们需要考虑的是我们通常用的drand48()是返回一个均匀分布的随机值，如何获得想要的分布的随机值呢？

我们先从正态分布的例子：

假设我们有一组数据，火者到站的时间与预估的到站时间差符合正态分布。我们这么利用drand48()获得正态分布的随机值呢？我们可以想到正态分布的CDF函数，因为其单调递增且值域为[0,1]：

我们的均匀随机分布正好可以是其CDF的y值，我们用y取计算X，就可以获得我们需要的正态分布(可以是任意)的随机值了。利用y计算X有两种方法，一种是直接利用公式反推的方法，另一种首先用程序根据pdf计算CDF表(离散形式)，然后利用drand48函数获取随机值y并在CDF表中索引并计算出X。示例代码如下：

#include <cstdlib>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm> // 正态分布pdf
float pdf(const float &x)
{ return 1 / sqrtf(2 * M_PI) * exp(-x * x * 0.5); } float sample(float *cdf, const uint32_t &nbins, const float &minBound, const float &maxBound)
{ float r = drand48(); float *ptr = std::lower_bound(cdf, cdf + nbins + 1, r); int off = std::max(0, (int)(ptr - cdf - 1)); float t = (r - cdf[off]) / (cdf[off + 1] - cdf[off]); float x = (off + t) / (float)(nbins); return minBound + (maxBound - minBound) * x;
} int main(int argc, char ** argv)
{ srand48(13); // 生成CDF表int nbins = 32; float minBound = -5, maxBound = 5; float cdf[nbins + 1], dx = (maxBound - minBound) / nbins, sum = 0; cdf[0] = 0.f; for (int n = 1; n < nbins; ++n) { float x = minBound + (maxBound - minBound) * (n / (float)(nbins)); float pdf_x = pdf(x) * dx; cdf[n] = cdf[n - 1] + pdf_x; sum += pdf_x; } cdf[nbins] = 1; // 我们的模拟int numSims = 100000; int numBins = 100; // to collect data on our sim int bins[numBins]; // to collect data on our sim memset(bins, 0x0, sizeof(int) * numBins); // 设定初始值0 const float dist = 10; // 10 km for (int i = 0; i < numSims; ++i) { float diff = sample(cdf, nbins, minBound, maxBound); // 利用y获取Xint whichBin = (int)(numBins * (diff - minBound) / (maxBound - minBound)); bins[whichBin]++; float time = 30 + diff; float speed = 60 * dist / time; } sum = 0; for (int i = 0; i < numBins; ++i) { float r = bins[i] / (float)numSims; printf("%f %f\n", 5 * (2 * (i /(float)(numBins)) -1), r); sum += r; } fprintf(stderr, "sum %f\n", sum); return 0;
}

本文不采用上述代码中的方法(上述代码的方法只是作为一个实例)，而是直接使用公式反推。我们的公式：

所以我们可以用drand48()获取 $\xi$ 的值，然后计算x，这样x就是满足pdf为 ${p}'(x)=8x/\pi^{2}$ 分布的随机值了。我们把使用了重要性采样和使用均匀分布的积分作对比：

#include <cstdio>
#include <cstdlib>
#include <cmath>
int main(int argc, char **argv)
{ srand48(13); int N = 16; for (int n = 0; n < 10; ++n) { float sumUniform = 0, sumImportance = 0; for (int i = 0; i < N; ++i) { float r = drand48(); sumUniform += sin(r * M_PI * 0.5); float xi = sqrtf(r) * M_PI * 0.5; // this is our X_i sumImportance += sin(xi) / ((8 * xi) / (M_PI * M_PI)); } sumUniform *= (M_PI * 0.5) / N; sumImportance *= 1.f / N; printf("%f %f\n", sumUniform, sumImportance); } return 0;
}

我们输出的结果如下表：

重要性采样的方差有明显的减少。

重要性采样(Importance Sampling)简介和简单样例实现相关推荐

重要性采样(Importance Sampling)详细学习笔记
重要性采样(Importance Sampling)详细学习笔记文章目录重要性采样(Importance Sampling)详细学习笔记前言: 参考主体: on-policy 和 off-pol ...
图形学数学基础之重要性采样(Importance Sampling)
作者:i_dovelemon 日期:2017/08/06 来源:CSDN 主题:Importance Sampling, PDF, Monte Carlo 引言前面的文章[图形学数学基础之基本蒙特卡 ...
matlab重要性采样,Importance Sampling (重要性采样)介绍 | 文艺数学君
摘要这一篇是关于重要性抽样(importance sampling)的介绍, 包括他的背景知识, 相关的数学转换和最后的例子. 简介重要性抽样(importance sampling)是一种近似的抽 ...
重要性采样Importance Sampling
参考: https://zhuanlan.zhihu.com/p/41217212 https://zhuanlan.zhihu.com/p/78720910?utm_source=wechat_se ...
重要性采样(importance sampling)
重要性采样是统计学习中一种常用的方法.在强化学习中通常和蒙特卡洛方法结合使用. 重要性采样是,使用另外一种分布来逼近所求分布一种方法. 具体形式是这样的:假设我们在想要求取目标分布PPP下函数f(x) ...
JDBC 连接Hive 简单样例（开启Kerberos）
今天在移动的云平台上通过jdbc连接hive,发现云平台使用了 kerberos的认证.与宁波实验环境不同. 发现一文解决了问题,转载如下: 原文地址:http://blog.csdn.net/zen ...
C语言单元测试之安装gtest教程及一个简单样例
准备工作安装包:gtest1.7.0版本(最新的1.8.0版本一直安装失败,1.7.0版本一次成功) 安装链接:百度网盘 https://pan.baidu.com/s/1mDy9sB3sBIMei ...
K8S Yaml 详细说明及简单样例
一.K8S Yaml 配置文件主要分为基本标签.元数据标签.资源内容 3 个部分基本标签 apiVersion: v1 #必选,版本号,例如v1 kind: Pod #必选,Pod 元数据标签 me ...
NASBench101-安装及简单样例使用指南
NASBench101-安装及简单样例使用指南 github地址:https://github.com/google-research/nasbench paper原文地址:https://arxiv ...

重要性采样(Importance Sampling)简介和简单样例实现

重要性采样(Importance Sampling)简介和简单样例实现

重要性采样(Importance Sampling)简介和简单样例实现相关推荐

最新文章

热门文章