c7-8 (779522)
Текст из файла
316Chapter 7.Random Numbers7.8 Adaptive and Recursive Monte CarloMethodsImportance SamplingThe use of importance sampling was already implicit in equations (7.6.6) and (7.6.7).We now return to it in a slightly more formal way. Suppose that an integrand f can be writtenas the product of a function h that is almost constant times another, positive, function g. Thenits integral over a multidimensional volume V isZZZf dV = (f /g) gdV = h gdV(7.8.1)In equation (7.6.7) we interpreted equation (7.8.1) as suggesting a change of variable toG, the indefinite integral of g.
That made gdV a perfect differential. We then proceededto use the basic theorem of Monte Carlo integration, equation (7.6.1). A more generalinterpretation of equation (7.8.1) is that we can integrate f by instead sampling h — not,however, with uniform probability density dV , but rather with nonuniform density gdV .
Inthis second interpretation, the first interpretation follows as the special case, where the meansof generating the nonuniform sampling of gdV is via the transformation method, using theindefinite integral G (see §7.2).More directly, one can go back and generalize the basic theorem (7.6.1) to the caseof nonuniform sampling: Suppose that points xi are chosen within the volume V with aprobability density p satisfyingZp dV = 1(7.8.2)The generalized fundamental theorem is that the integral of any function f is estimated, usingN sample points xi , . . .
, xN , bys ZZfhf 2 /p2 i − hf /pi2fI ≡ f dV =±pdV ≈(7.8.3)ppNwhere angle brackets denote arithmetic means over the N points, exactly as in equation(7.6.2). As in equation (7.6.1), the “plus-or-minus” term is a one standard deviation errorestimate. Notice that equation (7.6.1) is in fact the special case of equation (7.8.3), withp = constant = 1/V .What is the best choice for the sampling density p? Intuitively, we have alreadyseen that the idea is to make h = f /p as close to constant as possible. We can be morerigorous by focusing on the numerator inside the square root in equation (7.8.3), which isthe variance per sample point.
Both angle brackets are themselves Monte Carlo estimatorsof integrals, so we can write 2 2 Z 2Z2 Z 2Z2fffffS≡≈pdV−=(7.8.4)pdVdV−fdV−p2pp2ppWe now find the optimal p subject to the constraint equation (7.8.2) by the functional variation!Z2ZZ 2fδ(7.8.5)dV −f dV + λ p dV0=δppSample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited.
To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).This section discusses more advanced techniques of Monte Carlo integration. Asexamples of the use of these techniques, we include two rather different, fairly sophisticated,multidimensional Monte Carlo codes: vegas [1,2] , and miser [3]. The techniques that wediscuss all fall under the general rubric of reduction of variance (§7.6), but are otherwisequite distinct.7.8 Adaptive and Recursive Monte Carlo Methods317with λ a Lagrange multiplier.
Note that the middle term does not depend on p. The variation(which comes inside the integrals) gives 0 = −f 2 /p2 + λ or|f ||f |p= √ = R|f | dVλ(7.8.6)One curiosity is that one can add a constant to the integrand to make it all of one sign,since this changes the integral by a known amount, constant × V . Then, the optimal choiceof p always gives zero variance, that is, a perfectly accurate integral! The resolution ofthis seeming paradox (already mentioned at the endR of §7.6) is that perfect knowledge of pin equation (7.8.6) requires perfect knowledge of |f |dV , which is tantamount to alreadyknowing the integral you are trying to compute!If your function f takes on a known constant value in most of the volume V , it iscertainly a good idea to add a constant so as to make that value zero.
Having done that, theaccuracy attainable by importance sampling depends in practice not on how small equation(7.8.7) is, but rather on how small is equation (7.8.4) for an implementable p, likely only acrude approximation to the ideal.Stratified SamplingThe idea of stratified sampling is quite different from importance sampling. Let usexpand our notation slightly and let hhf ii denote the true average of the function f overthe volume V (namely the integral divided by V ), while hf i denotes as before the simplest(uniformly sampled) Monte Carlo estimator of that average:Z11 Xhhf ii ≡f dVhf i ≡f (xi)(7.8.8)VN iThe variance of the estimator, Var (hf i), which measures the square of the error of theMonte Carlo integration, is asymptotically related to the variance of the function, Var (f ) ≡hhf 2ii − hhf ii2, by the relationVar (hf i) =Var (f )N(7.8.9)(compare equation 7.6.1).Suppose we divide the volume V into two equal, disjoint subvolumes, denoted a and b,and sample N/2 points in each subvolume.
Then another estimator for hhf ii, different fromequation (7.8.8), which we denote hf i0 , is1(7.8.10)hf i0 ≡hf ia + hf ib2in other words, the mean of the sample averages in the two half-regions. The variance ofestimator (7.8.10) is given by 1Var hf i0 =Var hf ia + Var hf ib4Varb (f )1 Vara (f )(7.8.11)+=4N/2N/21=[Vara (f ) + Varb (f )]2NSample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use.
Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).where λ has been chosen to enforce the constraint (7.8.2).If f has one sign in the region of integration, then we get the obvious result that theoptimal choice of p — if one can figure out a practical way of effecting the sampling — isthat it be proportional to |f |.
Then the variance is reduced to zero. Not so obvious, but seento be true, is the fact that p ∝ |f | is optimal even if f takes on both signs. In that case thevariance per sample point (from equations 7.8.4 and 7.8.6) isZ2 Z2S = Soptimal =|f | dV−f dV(7.8.7)318Chapter 7.Random NumbersHere Vara (f ) denotes the variance of f in subregion a, that is, hhf 2iia − hhf ii2a , andcorrespondingly for b.From the definitions already given, it is not difficult to prove the relation11(7.8.12)[Vara (f ) + Varb (f )] + (hhf iia − hhf iib)224(In physics, this formula for combining second moments is the “parallel axis theorem.”)Comparing equations (7.8.9), (7.8.11), and (7.8.12), one sees that the stratified (into twosubvolumes) sampling gives a variance that is never larger than the simple Monte Carlo case— and smaller whenever the means of the stratified samples, hhf iia and hhf iib, are different.We have not yet exploited the possibility of sampling the two subvolumes with differentnumbers of points, say Na in subregion a and Nb ≡ N − Na in subregion b.
Let us do sonow. Then the variance of the estimator is1 Vara (f )Varb (f )Var hf i0 =+(7.8.13)4NaN − NaVar (f ) =Naσa=Nσa + σb(7.8.14)Here we have adopted the shorthand notation σa ≡ [Vara (f )]1/2 , and correspondingly for b.If Na satisfies equation (7.8.14), then equation (7.8.13) reduces to (σa + σb )2Var hf i0 =(7.8.15)4NEquation (7.8.15) reduces to equation (7.8.9) if Var (f ) = Vara (f ) = Varb (f ), in which casestratifying the sample makes no difference.A standard way to generalize the above result is to consider the volume V divided intomore than two equal subregions. One can readily obtain the result that the optimal allocation ofsample points among the regions is to have the number of points in each region j proportionalto σj (that is, the square root of the variance of the function f in that subregion). In spacesof high dimensionality (say d >∼ 4) this is not in practice very useful, however.
Dividing avolume into K segments along each dimension implies K d subvolumes, typically much toolarge a number when one contemplates estimating all the corresponding σj ’s.Mixed StrategiesImportance sampling and stratified sampling seem, at first sight, inconsistent with eachother. The former concentrates sample points where the magnitude of the integrand |f | islargest, that latter where the variance of f is largest. How can both be right?The answer is that (like so much else in life) it all depends on what you know and howwell you know it. Importance sampling depends on already knowing some approximation toyour integral, so that you are able to generate random points xi with the desired probabilitydensity p.
Характеристики
Тип файла PDF
PDF-формат наиболее широко используется для просмотра любого типа файлов на любом устройстве. В него можно сохранить документ, таблицы, презентацию, текст, чертежи, вычисления, графики и всё остальное, что можно показать на экране любого устройства. Именно его лучше всего использовать для печати.
Например, если Вам нужно распечатать чертёж из автокада, Вы сохраните чертёж на флешку, но будет ли автокад в пункте печати? А если будет, то нужная версия с нужными библиотеками? Именно для этого и нужен формат PDF - в нём точно будет показано верно вне зависимости от того, в какой программе создали PDF-файл и есть ли нужная программа для его просмотра.















