首页 > 代码库 > 以均匀随机分布生成任意斜率随机线性分布
以均匀随机分布生成任意斜率随机线性分布
本文的目的是利用(0,1)上的均匀分布随机数生成器生成区间为(imin,imax),斜率为slope的任意归一化线性随机数生成器。
借助(0,1)上的均匀随机数生成器,可以通过反函数法生成任意分布的随机数生成器。
对于C++,生成(imin,imax)上随机数生成器的代码为
double uniform_dist(double imin,double imax) { int temp; while((temp=rand())==RAND_MAX){ ; } return temp/RAND_MAX*(imax-imin)+imin); }
其中imin为随机变量最小值,imax为随机变量最大值。
在调用时,只需
x=uniform_dist(0,1);
就可以生成$(0,1)$上的任意随机变量。
对任累积意分布函数$F$,其概率密度为$f(x)$,则有
$F(x_1)=P(x\le x_1)=\int _{-\infty} ^{x_1} f(x‘) dx$
由累积分布函数的值范围$F \in (0,1),则若$F$为$(0,1)$上的随机变量,即可通过$F(x_1)$的反函数来求得与此对应的随机变量$x_1$。
因此要生成概率密度函数为$f(x)=slope*x+b$区间为$(imin,imax)$的随机变量,分为2步
1. 通过归一化求出b的值
$ \int _{-\infty} ^{\infty} f(x)dx=\int _{imin} ^{imax} f(x)dx=1$
可以求出
$b=\frac {1-\frac {1}{2} slope (imax^2-imin^2)}{imax-imin}$
2. 由累积分布函数
$F(x_1)= \int _{imin} ^{x_1} f(x)dx=F_1(x_1)-F_1(imin)$
其中$F_1$为概率密度函数$f(x)$的原函数。由$F(x_1)$的反函数$F^{-1} (x_1)$可以求得的随机变量$x_1 \in (imin,imax)$,并且服从概率密度函数$f(x)$。
带入方程得
$\frac {1}{2}a x_1^2+\frac {1-\frac {1}{2} slope (imax^2-imin^2)}{imax-imin}x_1 \\ -\frac {1}{2}\ slope\ imin^2-\frac {1-\frac {1}{2} slope (imax^2-imin^2)}{imax-imin} imin - F(x_1)=0$
令
$ a=\frac {1}{2}slope$
$ b=\frac {1-\frac {1}{2}\ slope \ (imax^2-imin^2)}{imax-imin}$
$ c=-\frac {1}{2}\ slope\ imin^2-\frac {1-\frac {1}{2} slope (imax^2-imin^2)}{imax-imin} imin - F(x_1)$
解得
$ x_1=\frac {-b + \sqrt {b^2-4ac}}{2a}$
C++实现程序如下
#include<iostream> #include<vector> #include<time.h> #include<math.h> #include<stdlib.h> using namespace std;
//*******************************************************************
// uniform distribution
//*******************************************************************
double uniform_dist(double imin,double imax) { int temp; while((temp=rand())==RAND_MAX){ ; } return((double) temp/RAND_MAX*(imax-imin)+imin); }
//*******************************************************************
// linear distribution for any slope
//*******************************************************************
double linear_dist(double a,double b,double imin,double imax,double slope) { double c; c=-0.5*slope*imin*imin-b*imin-uniform_dist(0,1); return ((-b+sqrt(b*b-4*a*c))/(2a)); }
int main() {
//===========================================================
// sequence to be sampled
//===========================================================
const N=1000000; vector< double> x(N); srand((unsigned)time(NULL));
//===========================================================
// parameters of the linear distribution
//===========================================================
double imin,imax,slope; imin=0.0; imax=3.0; slope=-0.2; double a,b,criterion; a=0.5*slope; b=(1-0.5*slope*(imax*imax-imin*imin))/(imax-imin);
//*********************************************************
// judge the reasonability of the linear distribution
//*********************************************************
if(slope>0){ criterion=slope*imin+b; if(criterion<0){ cerr<<"The probability of the variable must be larger than 0"<<endl; exit(1); } } else if(slope<0){ criterion=slope*imax+b; if(criterion<0){ cerr<<"The probability of the variable must be larger than 0"<<endl; exit(1); } }
//**********************************************************
// sampling the random data
//**********************************************************
for(int i=0;i!=x.size();++i){ x[i]=linear_dist(a,b,imin,imax,slope); }
//**********************************************************
// output the data
//**********************************************************
for(int i=0;i!=x.size();++i){ cout<<x[i]<<""; } return 0;
注意:
对斜率的选取和随机变量区间的选取需要考虑归一化的限制,也即是需要考虑每一个随机变量的概率要大于0。
以均匀随机分布生成任意斜率随机线性分布