重要抽样法（lmportanceSampling）●公式表示为：(n)da=4/)=含77J0，重要抽样法无疑是蒙特卡洛计算中最基本和常用的技巧之一。它在提高计算速度和增加数值结果的稳定性方面都有很大的潜力。它的局限性是能寻找出某分布密度函数g（x），并能解析求出其对应的分布函数G（x）的情况并不多。。例如，计算积分：exp(-/2)da该式当然可以解析算出，这单主要举例说明重要抽样法的求解过程式中的指数项随x增加而迅速减小，g（x）必须是另一有相同行为的函数，我们用指数函数在x=0附近的渐近展开来近似：
重要抽样法 (Importance Sampling)
重要抽样法2223crg(a2X22848由反函数法抽样：(1-a/2)da=c-2/4Glrn=2(1±V1-$)由于G(O）=0，故上式中括号内应取负号，即：n= 2(1 - V1 - s)原积分即可表示成：N(1- V1-)expV1-sii=1
重要抽样法
#include<iostream>// int 0^1 exp(-x/2) dx =2(1-exp(-0.5))#include<cstdlib>// importance sampling:g(x)=1-x/2.0#include<ctime>#include<iostream>#include<cmath>#include<cstdlib>#include<ctime>usingnamespacestddouble func(double x);#include<cmath>int main()usingnamespacestd;[srand(unsigned(time(O)));//seedtherandomizerdouble func(double x);inti, n;int main()doublecrude_mc,x,sum_sigma,fx,variance;{srand(unsigned(time(O));//seedtherandomizerinti, n;cin >> n;doublecrude_mc,tx,x,sum_sigma,fx,variance;crude mc=sum sigma=0.;//acrudeMonte-Carlomethodcin>> n;for(i= 1;i<= n; i++)crudemc=sumsigma=0.;x=rand()/(double)RANDMAXfor(i= 1;i<= n; i++)(fx=func(x);tx=3.0/4.0*rand()/(double)RANDMAX；crude mc+=fx//0=<tx<=3/4;sum_sigma+=fx*fx;x=2.*(1-sqrt(1.0-tx));Read inthe number of Monte-Carlo samples1000crude_mc=crude_mc/((double)n);fx=3.0/4.0*func(x);variance=0.00178519Integral 0.788191Exact=e.786939crude mc+=fx;sum_sigma=sum_sigma/((double)n);variance=sum_sigma-crude_mc*crude_mc;sum_sigma+=fx*fx;}variance=sqrt(variance)/double(n);)crude_mc= crude_mc/(double)n);sumsigma=sumsigma/((double)n);double func(double xRead in the number of Monte-Carlo samplesvariance=sumsigma-crude_mc*crude_mc;10001variance= o.0129141 Integral =0.782019variance=sqrt(variance)/double(n);)Exact=@.786939doublevaluedouble func(double x)value = exp(-x/2.0);( doublevalue;value = exp(-x/2.0)/(1.-x/2.0);return value;l//endoffunctiontoevaluatereturn value;J//endoffunctiontoevaluate
#include <iostream> #include <cstdlib> #include <ctime> #include <cmath> using namespace std; double func(double x); int main() { srand(unsigned(time(0))); // seed the randomizer int i, n; double crude_mc, x, sum_sigma, fx, variance; cin >> n; crude_mc = sum_sigma=0. ; // a crude Monte-Carlo method for ( i = 1; i <= n; i++){ x=rand()/(double)RAND_MAX; fx=func(x); crude_mc += fx; sum_sigma += fx*fx; } crude_mc = crude_mc/((double) n ); sum_sigma = sum_sigma/((double) n ); variance=sum_sigma-crude_mc*crude_mc; variance=sqrt(variance)/double(n); } double func(double x) { double value; value = exp(-x/2.0); return value; } // end of function to evaluate // int_0^1 exp(-x/2) dx =2(1-exp(-0.5)) // importance sampling: g(x)=1-x/2.0 #include <iostream> #include <cstdlib> #include <ctime> #include <cmath> using namespace std; double func(double x); int main() { srand(unsigned(time(0))); // seed the randomizer int i, n; double crude_mc, tx, x, sum_sigma, fx, variance; cin >> n; crude_mc = sum_sigma=0. ; for ( i = 1; i <= n; i++){ tx=3.0/4.0*rand()/(double)RAND_MAX; //0=<tx<=3/4; x=2.*(1-sqrt(1.0-tx)); fx=3.0/4.0*func(x); crude_mc += fx; sum_sigma += fx*fx; } crude_mc = crude_mc/((double) n ); sum_sigma = sum_sigma/((double) n ); variance=sum_sigma-crude_mc*crude_mc; variance=sqrt(variance)/double(n);} double func(double x) { double value; value = exp(-x/2.0)/(1.-x/2.0); return value; } // end of function to evaluate
//Breit-WignerSk-m?(sk-m) +mTede,=tan//mk=10.,gk=1.0;ds.mkmIk// int_0^200 sqrt(x)/(x-mk*mk)*(x-mk*mk)+mk*mk*gk*gk);// theta:atan(-10.)to atan(10.)intmain()int main()(srand(unsigned(time(0)); // seedtherandomizerfsrand(unsigned(time(O)));// seedtherandomizerinti, n;int i, n;doublecrude_mc,x,sum_sigma,fx,variance;doublecrude_mc,x,sum_sigma,fx,variance;cin>>n;cin >>n;crude mc=sum sigma=0.;doublepi=3.1415926535897crude mc=sum sigma=0.;for(i=1;i<=n;i++)for(i=1;i<=n;i++)x=rand()/(double)RANDMAXx=rand()/(double)RAND MAX;x=200.*x; //0-200x=-atan(10.)+2.*atan(10.)*x; // atan(-10)to atan(10)fx=func(x);fx=func(x);crude mc +=200.*fx;crude mc += 2.*atan(10.)*fxsumsigma+=200.*200.*fx*fx;!sum_sigma+=2.*atan(10.)*2.*atan(10.)*fx*fx;)crudemc=crudemc/((double)n);crude_mc = crude_mc/(double) n);sum_sigma=sum_sigma/(double)n);sum_sigma= sum_sigma/((double) n);variance=sum_sigma-crude_mc*crude_mcvariance=sum_sigma-crude_mc*crude_mc;variance=sqrt(variance)/double(n);variance=sqrt(variance)/double(n);cout<<"variance="<<variance<<" Integral ="cout<<"variance="<<variance<<" Integral ="<<crude_mc<<endl;<<crude_mc<<endl;1//endofmainprogram1//endofmainprogramdouble func(double x)//thisfunctiondefinesthefunctiontointegrate1double func(double x)(double value;double value;double mk=10., gk=1.0;doublemk=10..gk=1.0doubleax=tan(x)*mk*gk+mk*mk;value=sqrt(x)/((x-mk*mk)*(x-mk*mk)+mk*mk*gk*gk);value = sqrt(ax)/mk/gk;return value;)returnvalue:l//endoffunctiontoevaluate
int main() { srand(unsigned(time(0))); // seed the randomizer int i, n; double crude_mc, x, sum_sigma, fx, variance; cin >> n; crude_mc = sum_sigma=0. ; for ( i = 1; i <= n; i++){ x=rand()/(double)RAND_MAX; x=200.*x; // 0-200 fx=func(x); crude_mc += 200.*fx; sum_sigma += 200.*200.*fx*fx; } crude_mc = crude_mc/((double) n ); sum_sigma = sum_sigma/((double) n ); variance=sum_sigma-crude_mc*crude_mc; variance=sqrt(variance)/double(n); cout << " variance= " << variance << " Integral = " << crude_mc<< endl; } // end of main program // this function defines the function to integrate double func(double x) { double value; double mk=10., gk=1.0; value = sqrt(x)/((x-mk*mk)*(x-mk*mk)+mk*mk*gk*gk); return value;} int main() { srand(unsigned(time(0))); // seed the randomizer int i, n; double crude_mc, x, sum_sigma, fx, variance; cin >> n; crude_mc = sum_sigma=0. ; double pi=3.1415926535897; for ( i = 1; i <= n; i++){ x=rand()/(double)RAND_MAX; x=-atan(10.)+2.*atan(10.)*x; // atan(-10) to atan(10) fx=func(x); crude_mc += 2.*atan(10.)*fx; sum_sigma += 2.*atan(10.)*2.*atan(10.)*fx*fx; } crude_mc = crude_mc/((double) n ); sum_sigma = sum_sigma/((double) n ); variance=sum_sigma-crude_mc*crude_mc; variance=sqrt(variance)/double(n); cout << " variance= " << variance << " Integral = " << crude_mc<< endl; } // end of main program double func(double x) { double value; double mk=10., gk=1.0; double ax=tan(x)*mk*gk+mk*mk; value = sqrt(ax)/mk/gk; return value; } // end of function to evaluate // Breit-Wigner // mk=10., gk=1.0; // int_0^200 sqrt(x)/((x-mk*mk)*(x-mk*mk)+mk*mk*gk*gk); // theta: atan(-10.) to atan(10.)
The idea of stratified sampling is quite differentfrom importanceStratified抽样法sampling,bydividingvolumeVinto sub-regions.The optimalallocationistohavethenumberof pointsineachregionjproportional to gjVar(f)()=()d《J》三Var((f)) =NSuppose wedivide thevolumeVintotwo equal,disjoint subvolumes,denoted aand band sampleN/2pointsin each subvolume.Then another estimatorfor《f》,differentfromequation<f)",whichwedenote《f),is<fY'=(f)a+(f).)in other words,the mean of the sample averages in the two half-regions.The variance ofestimatoris givenby<f)"Var (<f)) =[Var ((f)a) + Var (《f>,)][Vara (f)Var, (f)14N/2N/21[Vara (f)+Var (f)]2NStratified后的方差总From the definitions already given, it is not difficult to prove the relation是小于之前的Var(f) = [Vara (f) + Var (f) + 1 (《f》>a - 《f)b)2
Stratified抽样法 The idea of stratified sampling is quite different from importance sampling, by dividing volume V into sub-regions. The optimal allocation is to have the number of points in each region j proportional to 𝜎𝑗 Stratified后的方差总 是小于之前的