随机采样之接受拒绝采样

news2024/11/25 0:45:18

之前提到的逆变换采样(Inverse Transform Sampling)是一种生成随机样本的方法,它利用累积分布函数(CDF)的逆函数来生成具有特定分布的随机变量。以下是逆变换采样的缺点:

  1. 计算复杂性:对于某些分布,找到累积分布函数(CDF)的逆函数可能是困难的,甚至是不可能的。
  2. 效率问题:对于具有重尾分布的随机变量,逆变换采样可能非常低效,因为CDF的逆可能需要大量的计算。
  3. 数值稳定性:在数值计算中,由于浮点数的精度限制,逆变换采样可能会引入误差,尤其是在CDF的值接近1时。

一、接受拒绝采样

接受-拒绝采样(Accept-Reject Sampling)方法是一种更为通用的采样方法,它可以用来生成具有任意分布的随机样本。这种方法不要求我们知道CDF的逆,而是利用一个简单的概率分布(称为提议分布)来生成样本,然后以一定的概率接受或拒绝这些样本。

接收-拒绝采样的基本步骤:

  1. 选择提议分布 g ( x ) g(x) g(x):选择一个容易从中抽样的分布 g ( x ) g(x) g(x),并且确保对于所有的 x x x,有 f ( x ) ≤ M ⋅ g ( x ) f(x) \leq M \cdot g(x) f(x)Mg(x),其中 f ( x ) f(x) f(x)是目标分布, M M M是一个正常数。

  2. 抽样:从提议分布 g ( x ) g(x) g(x)中抽取样本 x x x和从均匀分布 U ( 0 , 1 ) U(0, 1) U(0,1)中抽取样本 u u u

  3. 接受-拒绝条件:如果 u ≤ f ( x ) M ⋅ g ( x ) u \leq \frac{f(x)}{M \cdot g(x)} uMg(x)f(x),则接受 x x x作为目标分布 f ( x ) f(x) f(x)的一个样本;否则拒绝 x x x

接受拒绝采样可以使用下图进行表示(图片来源:【数之道】马尔可夫链蒙特卡洛方法是什么?十五分钟理解这个数据科学难点)。
在这里插入图片描述

二、接受拒绝采样证明

要证明接收-拒绝采样确实产生服从目标分布 f ( x ) f(x) f(x)的样本,我们需要证明对于所有的 x x x,有:
P ( X = x ) = f ( x ) (1) P(X=x) = f(x)\tag1 P(X=x)=f(x)(1)

其中 P ( X = x ) P(X=x) P(X=x)是样本 x x x被接受的概率。

证明:
  1. 接受概率:样本 x x x被接受的概率是 f ( x ) M ⋅ g ( x ) \frac{f(x)}{M \cdot g(x)} Mg(x)f(x),因为 u u u是从 U ( 0 , 1 ) U(0, 1) U(0,1)中抽取的。

  2. 联合概率:样本 x x x从提议分布 g ( x ) g(x) g(x)中抽取的概率是 g ( x ) g(x) g(x),并且 u u u [ 0 , f ( x ) M ⋅ g ( x ) ) [0, \frac{f(x)}{M \cdot g(x)}) [0,Mg(x)f(x))区间的概率是 f ( x ) M ⋅ g ( x ) \frac{f(x)}{M \cdot g(x)} Mg(x)f(x)。因此,联合概率是:

    P ( X = x , U ≤ f ( x ) M ⋅ g ( x ) ) = g ( x ) ⋅ f ( x ) M ⋅ g ( x ) = f ( x ) M (2) P(X=x, U \leq \frac{f(x)}{M \cdot g(x)}) = g(x) \cdot \frac{f(x)}{M \cdot g(x)} = \frac{f(x)}{M}\tag2 P(X=x,UMg(x)f(x))=g(x)Mg(x)f(x)=Mf(x)(2)

  3. 边缘概率:现在我们需要计算 X X X的边缘概率 P ( X = x ) P(X=x) P(X=x),即样本 x x x被接受的总概率。由于 u u u是均匀分布的,我们可以将联合概率在 u u u的所有可能值上积分:

    P ( X = x ) = ∫ 0 1 P ( X = x , U = u )   d u = ∫ 0 1 f ( x ) M   d u = f ( x ) M ⋅ ∫ 0 1 d u = f ( x ) M (3) P(X=x) = \int_0^1 P(X=x, U=u) \, du = \int_0^1 \frac{f(x)}{M} \, du = \frac{f(x)}{M} \cdot \int_0^1 du = \frac{f(x)}{M}\tag3 P(X=x)=01P(X=x,U=u)du=01Mf(x)du=Mf(x)01du=Mf(x)(3)

  4. 归一化常数:由于 M M M是使得 f ( x ) ≤ M ⋅ g ( x ) f(x) \leq M \cdot g(x) f(x)Mg(x)对所有 x x x成立的最小常数,我们可以将上式中的 M M M移到 f ( x ) f(x) f(x)的定义中,从而得到:

    P ( X = x ) = f ( x ) (4) P(X=x) = f(x)\tag4 P(X=x)=f(x)(4)
    这就证明了接收-拒绝采样确实产生了服从目标分布 f ( x ) f(x) f(x)的样本。

三、接受拒绝采样模拟

借用作者anshuai_aw1的例子,设我们需要采样的pdf为:
f ( x ) = 0.3 exp ⁡ ( − ( x − 0.3 ) 2 ) + 0.7 exp ⁡ ( − ( x − 2 ) 2 / 0.3 ) (5) f(x)=0.3 \exp \left(-(x-0.3)^{2}\right)+0.7 \exp \left(-(x-2)^{2} / 0.3\right)\tag5 f(x)=0.3exp((x0.3)2)+0.7exp((x2)2/0.3)(5)
其归一化常数为 Z = 1.2113 Z = 1.2113 Z=1.2113, 参考分布为 g ( x ) = N ( μ = 1.4 , σ 2 = ( 1. 2 2 ) ) g(x) =N(\mu=1.4,\sigma^2=(1.2^2)) g(x)=N(μ=1.4,σ2=(1.22)), M = 2.5 M=2.5 M=2.5, 以确保 M ⋅ g ( x ) ≥ f ( x ) M \cdot g(x) \geq f(x) Mg(x)f(x)。采样的代码如下:

import numpy as np
import matplotlib.pyplot as plt

def f(x):
    return (0.3*np.exp(-(x-0.3)**2) + 0.7* np.exp(-(x-2.)**2/0.3))/1.2113
x = np.arange(-4.,6.,0.01)
plt.plot(x,f(x),color = "red")

size = int(1e+07)
mu = 1.4
sigma = 1.2
M = 2.5

x = np.random.normal(loc = mu,scale = sigma, size = size)
g_x = 1/(np.sqrt(2*np.pi)*sigma)*np.exp(-0.5*(x-mu)**2/sigma**2)
u = np.random.uniform(low = 0, high = M*g_x, size = size)  #在[0,M*g_x]中均匀采样
fx =  0.3*np.exp(-(x-0.3)**2) + 0.7* np.exp(-(x-2.)**2/0.3)
sample = x[u <= fx] # u < fx(x)
plt.hist(sample,bins=150, density=True, edgecolor='black')
plt.show()

结果如下,其中红色曲线的是公式(5)所示pdf的图像,蓝色区域是采样结果,可见采样结果跟真实分布几乎一致。
在这里插入图片描述

参考资料:

[1]【数之道】马尔可夫链蒙特卡洛方法是什么?十五分钟理解这个数据科学难点
[2] 逆采样(Inverse Sampling)和拒绝采样(Reject Sampling)原理详解

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.coloradmin.cn/o/2238069.html

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈,一经查实,立即删除!

相关文章

软件设计师:排序算法总结

一、直接插入 排序方式&#xff1a;从第一个数开始&#xff0c;拿两个数比较&#xff0c;把后面一位跟前面的数比较&#xff0c;把较小的数放在前面一位 二、希尔 排序方式&#xff1a;按“增量序列&#xff08;步长&#xff09;”分组比较&#xff0c;组内元素比较交换 假设…

信息安全工程师(78)网络安全应急响应技术与常见工具

前言 网络安全应急响应是指为应对网络安全事件&#xff0c;相关人员或组织机构对网络安全事件进行监测、预警、分析、响应和恢复等工作。 一、网络安全应急响应技术 网络安全应急响应组织 构成&#xff1a;网络安全应急响应组织主要由应急领导组和应急技术支撑组构成。领导组负…

Kafka 的一些问题,夺命15连问

kafka-中的组成员 kafka四大核心 生产者API 允许应用程序发布记录流至一个或者多个kafka的主题&#xff08;topics&#xff09;。 消费者API 允许应用程序订阅一个或者多个主题&#xff0c;并处理这些主题接收到的记录流 StreamsAPI 允许应用程序充当流处理器&#xff08;s…

精选5款小程序设计工具,助力设计之路璀璨前行

在当今数字化浪潮中&#xff0c;小程序的重要性日益凸显&#xff0c;无论是电商、社交还是服务领域&#xff0c;小程序都成为连接用户与品牌的关键桥梁。而一款优秀的小程序离不开精心的设计&#xff0c;以下 5 款小程序设计工具将成为你设计事业的得力助手。 一、即时设计 即…

亚马逊评论爬虫+数据分析

爬取评论 做分析首先得有数据&#xff0c;数据是核心&#xff0c;而且要准确&#xff01; 1、爬虫必要步骤&#xff0c;选好框架 2、开发所需数据 3、最后测试流程 这里我所选框架是seleniumrequest&#xff0c;很多人觉得selenium慢&#xff0c;确实不快&#xff0c;仅针对此…

量子计算及其在密码学中的应用

&#x1f493; 博客主页&#xff1a;瑕疵的CSDN主页 &#x1f4dd; Gitee主页&#xff1a;瑕疵的gitee主页 ⏩ 文章专栏&#xff1a;《热点资讯》 量子计算及其在密码学中的应用 量子计算及其在密码学中的应用 量子计算及其在密码学中的应用 引言 量子计算概述 定义与原理 发展…

论文笔记:no pose,no problem-基于dust3r输出GS参数实现unpose稀疏重建

1.摘要 我们引入了 NoPoSplat&#xff0c;这是一种前馈模型&#xff0c;能够从未设置的稀疏多视图图像中重建由 3D 高斯参数化的 3D 场景。 我们的模型专门使用光度损失进行训练&#xff0c;在推理过程中实现了实时 3D 高斯重建。 为了消除重建过程中对准确pose的需要&#xff…

godot--自定义边框/选中时样式 StyleBoxTexture

前提知识&#xff1a; stylebox就像一个贴图&#xff0c;把图案贴到控件是。多个stylebox同时生效的话&#xff0c;那当然也有层级之分&#xff0c;上层覆盖下层&#xff08;可以设置透明度来显示下层&#xff09; 关于主题的概念&#xff1a; godot——主题、Theme、StyleB…

ReactPress 安装指南:从 MySQL 安装到项目启动

ReactPress Github项目地址&#xff1a;https://github.com/fecommunity/reactpress 欢迎Star。 ReactPress 是一个基于 React 的开源发布平台&#xff0c;适用于搭建博客、网站或内容管理系统&#xff08;CMS&#xff09;。本文将详细介绍如何安装 ReactPress&#xff0c;包括…

BMC运维管理:IPMI实现服务器远控制

IPMI实现服务器远控制 实操一、使用IPMI重置BMC用户密码实操二、使用IPMI配置BMC的静态IP实操三、IPMI实现BMC和主机控制操作实操四、ipmitool查看服务器基本信息实操五、ipmitool实现问题定位BMC(Baseboard Management Controller,基板管理控制器)是服务器硬件的一个独立管…

手机上用什么方法可以切换ip

手机上用什么方法可以切换IP&#xff1f;在某些特定情境下&#xff0c;用户可能需要切换手机的IP地址&#xff0c;以满足网络安全、隐私保护或绕过地域限制等需求。下面以华为手机为例&#xff0c;将详细介绍手机IP地址切换的几种方法&#xff0c;帮助用户轻松实现这一目标。 一…

鸿蒙多线程开发——并发模型对比(Actor与内存共享)

1、概 述 并发是指在同一时间段内&#xff0c;能够处理多个任务的能力。为了提升应用的响应速度与帧率&#xff0c;以及防止耗时任务对主线程的干扰&#xff0c;HarmonyOS系统提供了异步并发和多线程并发两种处理策略。 异步并发&#xff1a;指异步代码在执行到一定程度后会被…

【NLP】使用 PyTorch 从头构建自己的大型语言模型 (LLM)

读完这篇文章后&#xff0c;你会取得什么成就&#xff1f;你将能够自己构建和训练大型语言模型 (LLM)&#xff0c;同时与我一起编写代码。虽然我们正在构建一个将任何给定文本从英语翻译成马来语的 LLM&#xff0c;但你可以轻松地修改此 LLM 架构以用于其他语言翻译任务。 LLM…

css:还是语法

emmet的使用 emmet是一个插件&#xff0c;Emmet 是 Zen Coding 的升级版&#xff0c;由 Zen Coding 的原作者进行开发&#xff0c;可以快速的编写 HTML、CSS 以及实现其他的功能。很多文本编辑器都支持&#xff0c;我们只是学会使用它&#xff1a; 生成html结构 <!-- emme…

YOLO即插即用---PConv

Run, Don’t Walk: Chasing Higher FLOPS for Faster Neural Networks 论文地址&#xff1a; 1. 论文解决的问题 2. 解决问题的方法 3. PConv 的适用范围 4. PConv 在目标检测中的应用 5. 评估方法 6. 潜在挑战 7. 未来研究方向 8.即插即用代码 论文地址&#xff1a; …

小白NAS磁盘规划实践:一次科学、高效的存储旅程

引言 如今,数字化生活正逐步渗透到我们生活的方方面面。从家庭影音到工作文件,从珍贵的照片到大型游戏库,数据的存储需求日益增加。许多朋友开始关注NAS(网络附加存储)设备。作为一个专为数据存储和管理设计的系统,NAS能为我们提供安全、高效的存储方案。但如何科学地规…

ADC前端控制与处理模块--AD7606_Module

总体框架 AD7606_Module主要由3个模块组成组成&#xff0c;AD7606_Data_Pkt和AD7606_Drive以及AD7606_ctrl。 1.AD7606_Data_Pkt主要作用是把AD芯片数据组好数据包&#xff0c;然后发送给上位机&#xff1b; 2.AD7606_Drive主要负责和芯片的交互部分 3.AD7606_ctrl控制模块的作…

注意力机制篇 | YOLO11改进 | 即插即用的高效多尺度注意力模块EMA

前言&#xff1a;Hello大家好&#xff0c;我是小哥谈。与传统的注意力机制相比&#xff0c;多尺度注意力机制引入了多个尺度的注意力权重&#xff0c;让模型能够更好地理解和处理复杂数据。这种机制通过在不同尺度上捕捉输入数据的特征&#xff0c;让模型同时关注局部细节和全局…

dell服务器安装ESXI8

1.下载镜像在官网 2.打开ipmi&#xff08;idrac&#xff09;&#xff0c;将esxi镜像挂载&#xff0c;然后服务器开机 3.进入bios设置cpu虚拟化开启&#xff0c;进入boot设置启动选项为映像方式 4..进入安装引导界面3.加载完配置进入安装 系统提示点击继 5.选择安装磁盘进行…

Linux -- 操作系统(软件)

目录 什么是操作系统&#xff1f; 计算机的层状结构 为什么要有操作系统 操作系统到底层硬件 驱动程序 操作系统如何管理硬件&#xff1f; 操作系统到用户 系统调用接口 库函数 回到问题 什么是操作系统&#xff1f; 操作系统&#xff08;Operating System&#xf…