(一)连续随机量的生成-接受拒绝重采样

news/2024/2/29 3:15:05

连续随机量的生成-接受拒绝重采样

  • 1. 接受-拒绝重采样方法
  • 2. Python编程实现

1. 接受-拒绝重采样方法

重采样方法由两个步骤组成,第一个步骤提供近似分布的随机量,第二个步骤是校正机制。 在下文中,我们用 f f f 表示目标分布,用 g g g 表示辅助分布。

本课程我们只研究经典的接受-拒绝方法。

我们的目标是从概率密度为 f f f的分布中抽取一个随机量 X X X,这可能很难抽取。 我们引入辅助概率分布 g g g并采取以下假设:
(A) 假设存在一些 M > 0 M>0 M>0 使得
f ( x ) g ( x ) ≤ M 对于所有  x 。  \frac{f(x)}{g(x)} \leq M \quad \text { 对于所有 } x \text {。 } g(x)f(x)M 对于所有 x 

接受-拒绝方法的伪代码:

  • 从分布 g g g 中采样随机数量 Y
  • U ( [ 0 , 1 ] ) U([0,1]) U([0,1]) 中抽取随机数量 U U U
    • 如果 U ⩽ f ( Y ) M g ( Y ) U \leqslant \frac{f(Y)}{M g(Y)} UMg(Y)f(Y),接受,即。 X = Y X=Y X=Y
    • 如果 U > f ( Y ) M g ( Y ) U>\frac{f(Y)}{M g(Y)} U>Mg(Y)f(Y),则拒绝并重复该过程。

该算法的理论证明:

让我们考虑一维情况。 对于任何 x ∈ R x \in \mathbb{R} xR,考虑
P ( X ⩽ x ) = P ( Y ⩽ x ∣ f ( Y ) M g ( Y ) ⩾ U ) = P ( Y ⩽ x , f ( Y ) M g ( Y ) ⩾ U ) P ( f ( Y ) M g ( Y ) ⩽ U ) = ∫ − ∞ x ( ∫ 0 f ( Y ) M g ( y ) d u ) g ( y ) d y ∫ − ∞ + ∞ ( ∫ 0 f ( y ) M g ( y ) d u ) g ( y ) d y ( ∵ f r a c f ( y ) M g ( y ) ⩽ 1 对于所有  y ) = ∫ − ∞ x f ( y ) d y ∫ − ∞ + ∞ f ( y ) d y \begin{aligned} P(X\leqslant x) & =P\left(Y \leqslant x \mid \frac{f(Y)}{M g(Y)} \geqslant U\right) \\ & =\frac{P\left(Y \leqslant x, \frac{f(Y)}{M g(Y)} \geqslant U\right)}{P\left(\frac{f(Y)}{ M g(Y)} \leqslant U\right)} \\ & =\frac{\int_{-\infty}^x\left(\int_0^{\frac{f(Y)}{M g(y)}} d u\right) g(y) d y}{\int_ {-\infty}^{+\infty}\left(\int_0^{\frac{f(y)}{M g(y)}} d u\right) g(y) d y}\left(\because \ frac{f(y)}{M g(y)} \leqslant 1 \text { 对于所有 } y\right) \\ & =\frac{\int_{-\infty}^x f(y) d y}{\int_{-\infty}^{+\infty} f(y) d y} \end{aligned} P(Xx)=P(YxMg(Y)f(Y)U)=P(Mg(Y)f(Y)U)P(Yx,Mg(Y)f(Y)U)=+(0Mg(y)f(y)du)g(y)dyx(0Mg(y)f(Y)du)g(y)dy( fracf(y)Mg(y)1 对于所有 y)=+f(y)dyxf(y)dy
因此, X X X 具有密度为 f f f 的分布。

2. Python编程实现

接受-预测方法可用于对分布乘上一个常数的随机量进行采样。 一个示例是从具有以下密度的分布中采样随机量 X X X
1 C ⋅ 1 1 + ∣ x − 2 ∣ 3 \frac{1}{C} \cdot \frac{1}{1+|x-2|^3} C11+x231
其中 C = ∫ − ∞ + ∞ 1 1 + ∣ x − 2 ∣ 3 d x C=\int_{-\infty}^{+\infty} \frac{1}{1+|x-2|^3} d x C=+1+x231dx,无法显式计算出来。

编写伪代码,通过接受-拒绝方法从分布(1.3.1)中采样随机量(提示:取 M = 5 M=5 M=5 )。 编写一个程序来调整接受-拒绝过程 1000 次,计算创建了多少个随机量,并绘制直方图。

import numpy as np
import matplotlib.pyplot as plt# Define the target distribution function f(x)
def target_distribution(x):return 1 / (1 + np.abs(x - 2) ** 3)# Define the auxiliary distribution (Cauchy distribution)
def auxiliary_distribution(x, x0, gamma):return 1 / (np.pi * gamma * (1 + ((x - x0) / gamma) ** 2))# Set the number of samples to generate
num_samples = 1000# Initialize arrays to store accepted and rejected samples
accepted_samples = []
rejected_samples = []# Parameters for the Cauchy distribution
x0 = 2  # Center of the Cauchy distribution
gamma = 1  # Scale parameter of the Cauchy distribution# Set the maximum value for the acceptance ratio
max_f_x = 5# Acceptance rate
acceptance_rate = 0# Generate samples using accept-reject method
for _ in range(num_samples):# Generate a random sample from the Cauchy distributionx_sample = np.random.cauchy(x0, gamma)# Generate a uniform random number for acceptance/rejectionu = np.random.uniform(0, max_f_x)# Calculate the acceptance probabilityacceptance_prob = target_distribution(x_sample) / (auxiliary_distribution(x_sample, x0, gamma))# Accept or reject the sampleif u < acceptance_prob:accepted_samples.append(x_sample)acceptance_rate += 1else:rejected_samples.append(x_sample)# Calculate the acceptance rate
acceptance_rate /= num_samples# Plot the histogram of accepted samples
plt.hist(accepted_samples, bins=30, density=True, alpha=0.5, label='Accepted Samples')
plt.plot(x_range, target_distribution(x_range), 'r', label='Target Distribution')
plt.xlabel('x')
plt.ylabel('Density')
plt.legend()
plt.title(f'Acceptance Rate: {acceptance_rate:.2%}')
plt.show()print(f'Number of accepted samples: {len(accepted_samples)}')

http://www.ppmy.cn/news/1081316.html

相关文章

学信网学历电子注册备案表 下载方法

如果你需要成人高考的专升本 那么就会需要 学信网学历电子注册备案表 全称叫 教育部学历证书电子注册备案表 是学信网依托全国高等教育学生信息数据库&#xff0c;对学生的学历信息提供的在线验证报告&#xff0c;是我们验证学历真伪的一份报告 学历电子注册备案表在 考研、考…

2023-08-29力扣每日一题

链接&#xff1a; 823. 带因子的二叉树 题意&#xff1a; 用给的数字建二叉树&#xff0c;要求父节点是子节点的乘积 解&#xff1a; 乐了 1500ms30MB //注释版120ms18MB 实际代码&#xff1a; #include<bits/stdc.h> using namespace std; static constexpr int…

Linux知识点 -- Linux多线程(四)

Linux知识点 – Linux多线程&#xff08;四&#xff09; 文章目录 Linux知识点 -- Linux多线程&#xff08;四&#xff09;一、线程池1.概念2.实现3.单例模式的线程池 二、STL、智能指针和线程安全1.STL的容器是否是线程安全的2.智能指针是否是线程安全的 三、其他常见的各种锁…

安服面试 --- 01

1、常用渗透工具 burp、nmap、sqlmap、蚁剑、御剑、冰蝎、cobalt strike等 2、渗透测试中&#xff0c;拿到目标公司站点&#xff0c;接下来应该怎么做&#xff1f; &#xff08;1&#xff09;信息收集&#xff1a;收集目标公司的相关信息。包括域名、ip地址、子域名、开放端…

code 架构

目录 1. code 架构1.1. 代码质量的评判的维度1.2. 架构师1.3. 基础平台篇1. code 架构 1.1. 代码质量的评判的维度 可阅读性 (方便代码流转)可扩展性 / 可维护性(方便修改功能, 添加新功能)可测试性 (质量管理)可复用性 (简化后续功能开发的难度)1.2. 架构师 软件工程是一项非…

c#垃圾回收(Garbage Collection)

在C#中&#xff0c;垃圾回收&#xff08;Garbage Collection&#xff09;是一种自动管理内存的机制。它负责跟踪和释放不再使用的内存&#xff0c;以便程序可以有效地使用内存资源。 C#中的垃圾回收器是由.NET运行时&#xff08;CLR&#xff09;提供和管理的。它使用了一种叫做…

VB高校固定资产管理系统设计与实现

摘要 随着电脑的普及与使用,现在的管理也提升了一个档次,渐渐实现了无纸化办公,即从原来的人工记录管理模式转变为电脑一体化管理。高校是科研的阵地,后勤的高校固定资产管理系统也应该一改传统的人工管理,更加信息化,时代化,节省人力物力,提高效率。基于这一点,开发此…

Lesson3-5:OpenCV图像处理---模版匹配和霍夫变换

学习目标 掌握模板匹配的原理&#xff0c;能完成模板匹配的应用理解霍夫线变换的原理&#xff0c;了解霍夫圆检测知道使用OpenCV如何进行线和圆的检测 1 模板匹配 1.1 原理 所谓的模板匹配&#xff0c;就是在给定的图片中查找和模板最相似的区域&#xff0c;该算法的输入包括…

【Locomotor运动模块】抓取:按朝向抓取(Orientation Handler)案例

文章目录 案例原理 案例 左右手柄抓宝剑时&#xff0c;宝剑的朝向不同 L35 一个手柄对应一个抓取点 原理 1、左右手柄分别抓取的是宝剑上的不同抓取点——GenericOrientation Handle通用朝向把手 它是我们设置“按朝向抓取”&#xff08;Orientation Handler&#xff09;时&…

kvm虚拟机开启VNC功能

停止kvm虚拟机 virsh shutdown 虚拟机名称 在kvm虚拟机配置文件里面添加如下内容 <graphics typevnc port-1 autoportyes listen0.0.0.0 keymapen-us passwd123456> 启动kvm虚拟机 virsh start 虚拟机名称 得到虚拟机进程id ps -ef|grep 虚拟机名称 得到虚拟机vnc…

人工智能(AI)在材料科学方面的应用

人工智能&#xff08;AI&#xff09;在材料科学方面的应用日益增多&#xff0c;主要包括以下几个方面&#xff1a; 材料设计和发现&#xff1a;通过机器学习和深度学习算法&#xff0c;预测材料的性质和特性&#xff0c;在材料研究和开发中起到重要的作用。例如&#xff0c;使用…

固定资产制度怎么完善管理?

固定资产管理制度的完善管理可以从以下几个方面入手&#xff1a;  建立完善的资产管理制度&#xff0c;可以及时掌握企业资产的信息状况&#xff0c;使资产管理更加明确&#xff0c;防止资产流失。  加大固定资产监管力度&#xff0c;从配置资产、使用资产到处置资产进行全…

青翼科技基于VITA57.1的16路数据收发处理平台产品手册

FMC211是一款基于VITA57.1标准规范的实现16路LVDS数据采集、1路光纤数据收发处理FMC子卡模块。 该板卡支持2路CVBS&#xff08;复合视频&#xff09;视频输入&#xff0c;能够自动检测标准的模拟基带电视信号&#xff0c;并将其转变为8位ITU-R.656接口信号或者4:2:2分量视频信…

token无感刷新 前端实现 一看就会!

场景&#xff1a;在往常的接口中我们通过从后端获取的token来进行身份识别&#xff0c;并且获取一些登录后的才能获取的一些数据&#xff0c;那么这个token有效期很短的时候&#xff0c;那么总不能让用户再登录一次来刷新token&#xff0c;所以我们需要提供一项操作&#xff0c…

即时通讯开发应用中的实时消息推送技术

即时通讯开发领域正以前所未有的速度蓬勃发展&#xff0c;实时消息推送技术成为促进即时通讯应用体验的关键要素。本文将深入探讨即时通讯应用中的实时消息推送技术&#xff0c;为读者呈现这一领域的全貌。 2. 实时消息推送的重要性 在当今数字化时代&#xff0c;人们日益需要…

Java网络编程-Socket实现数据通信

文章目录 前言网络编程三要素IP地址和端口号传输协议Socket 使用Scoket实现网络通信TCPTCP通信-发送方TCP通信-接收方结果 UDPUDP通信-发送方UDP通信-接收方结果 总结 前言 本文主要是为下一篇Websockt做铺垫&#xff0c;大家了解socket的一些实现。 网络编程三要素 网络编程…

统一使用某一个包管理工具,比如yarn pnpm

原因&#xff1a;前端每个人的习性不一样&#xff0c;有人用npm 有人用yarn等包管理工具&#xff0c;混合下载插件容易出bug&#xff0c;就用个小工具锁住就行了&#xff0c;只能使用yarn或者pnpm反向下载依赖和下载插件。不然就报错 1.在项目主目录下创建preinstall.js // 如…

1.12 进程注入ShellCode套接字

在笔者前几篇文章中我们一直在探讨如何利用Metasploit这个渗透工具生成ShellCode以及如何将ShellCode注入到特定进程内&#xff0c;本章我们将自己实现一个正向ShellCodeShell&#xff0c;当进程被注入后&#xff0c;则我们可以通过利用NC等工具连接到被注入进程内&#xff0c;…

国标视频云服务EasyGBS国标视频平台迁移服务器后无法启动的问题解决方法

国标视频云服务EasyGBS支持设备/平台通过国标GB28181协议注册接入&#xff0c;并能实现视频的实时监控直播、录像、检索与回看、语音对讲、云存储、告警、平台级联等功能。平台部署简单、可拓展性强&#xff0c;支持将接入的视频流进行全终端、全平台分发&#xff0c;分发的视频…

Kotlin return 和 loop jump

再聊 return 在上一篇文章《Kotlin inline、noinline、crossinline 深入解析》 我们介绍到,在 lambda 中不能使用 return,除非该函数是 inline 的。如果该高阶函数是 inline ,调用该函数时,在传入的 lambda 中使用 return,则 return 的是离它最近的 enclosing function,…
最新文章