Marsaglia polära metod

Marsaglias normala polära metod är en pseudo-slumpmässigt talsamplingsmetod för att generera ett par oberoende standard slumpvariabler .

Standard normala slumpvariabler används ofta inom datavetenskap , beräkningsstatistik och i synnerhet i tillämpningar av Monte Carlo-metoden .

Den polära metoden fungerar genom att välja slumpmässiga punkter ( x , y ) i kvadraten −1 < x < 1, −1 < y < 1 tills

$0<s=x^{2}+y^{2}<1,\,$

och sedan returnera det nödvändiga paret av normala slumpvariabler som

$x{\sqrt {\frac {-2\ln(s)}{s}}}\,,\ \ y {\sqrt {\frac {-2\ln(s)}{s}}},$

eller på motsvarande sätt

${\frac {x}{\sqrt {s}}}{\sqrt {-2\ln(s)}} \,,\ \ {\frac {y}{\sqrt {s}}}{\sqrt {-2\ln(s)}},$

där $x/{\sqrt {s}}$ och $y/{\sqrt {s}}$ representerar cosinus och sinus för vinkeln som vektorn ( x , y ) gör med x -axeln.

Teoretisk grund

Den bakomliggande teorin kan sammanfattas enligt följande:

Om u är likformigt fördelad i intervallet 0 ≤ u < 1, så är punkten (cos(2π u ), sin(2π u )) likformigt fördelad på enhetens omkrets x ² + y ² = 1, och multiplicera den punkten med en oberoende stokastisk variabel ρ vars fördelning är

$\Pr(\rho <a)=\int _{0}^{a}re^{-r^{2}/ 2}\,dr$

kommer att producera en poäng

$\left(\rho \cos(2\pi u),\rho \sin(2\pi u)\right)$

vars koordinater är gemensamt fördelade som två oberoende standard normala slumpvariabler.

Historia

Denna idé går tillbaka till Laplace , som Gauss tillskriver att han hittat ovanstående

$I=\int _{-\infty }^{\infty }e^{-x^{2}/2}\,dx$

genom att ta kvadratroten av

$I^{2}=\int _{-\infty }^{\infty }\int _{-\infty }^{\infty }e^{-(x^{2}+y^{2 })/2}\,dx\,dy=\int _{0}^{2\pi }\int _{0}^{\infty }re^{-r^{2}/2}\,dr \,d\theta .$

Transformationen till polära koordinater visar att θ är likformigt fördelad (konstant densitet) från 0 till 2π, och att det radiella avståndet r har densitet

$re^{-r^{2}/2}.\,$

( r ² har lämplig chi-kvadratfördelning .)

Denna metod för att producera ett par oberoende standardnormalvariationer genom att radiellt projicera en slumpmässig punkt på enhetens omkrets till ett avstånd som ges av kvadratroten av en chi-kvadrat-2-variat kallas den polära metoden för att generera ett par normala slumpvariabler ,

Praktiska överväganden

En direkt tillämpning av denna idé,

$x={\sqrt {-2\ln (u_{1})}}\cos(2\pi u_{2}),\quad y={\sqrt {-2\ln(u_{1})}}\sin(2\pi u_{2} )$

kallas Box–Muller-transformen , där chi-variaten vanligtvis genereras som

${\sqrt {-2\ln(u_{1})}},$

men den transformationen kräver logaritm, kvadratrot, sinus och cosinusfunktioner. På vissa processorer kan cosinus och sinus för samma argument beräknas parallellt med en enda instruktion. Särskilt för Intel-baserade maskiner kan man använda fsincos assembler-instruktion eller expi-instruktionen (tillgänglig t.ex. i D), för att beräkna komplexa

$\operatorname {expi} (z)=e^{iz}=\cos(z)+i\sin (z),\,$

och bara separera de verkliga och imaginära delarna.

Obs: För att explicit beräkna den komplexa polära formen använd följande ersättningar i den allmänna formen,

Låt $r={\sqrt {-2\ln(u_{1})}}$ och $z=2\pi u_{2}.$ Sedan

$\ re^{iz}={\sqrt {-2\ln(u_{1})}}e^{i2\pi u_{2}}={\sqrt {-2\ln(u_{1 })}}\vänster[\cos(2\pi u_{2})+i\sin(2\pi u_{2})\höger].$

Däremot tar den polära metoden här bort behovet av att beräkna en cosinus och sinus. Istället, genom att lösa en punkt på enhetscirkeln, kan dessa två funktioner ersättas med x- och y -koordinaterna normaliserade till ${\sqrt {x^{2}+y^{2} }}$ radie. I synnerhet projiceras en slumpmässig punkt ( x , y ) inuti enhetscirkeln på enhetens omkrets genom att ställa in $s=x^{2}+y^{2}$ och bilda punkt

$\left({\frac {x}{\sqrt {s}}},{\frac {y}{\sqrt {s}}}\right),\,$

vilket är en snabbare procedur än att beräkna cosinus och sinus. Vissa forskare hävdar att den villkorade if-instruktionen (för att avvisa en punkt utanför enhetscirkeln), kan göra program långsammare på moderna processorer utrustade med pipelining och förgreningsprediktion. Även denna procedur kräver cirka 27 % fler utvärderingar av den underliggande slumptalsgeneratorn (endast $\pi /4\approx 79\%$ av genererade punkter ligger inuti enhetscirkeln).

Den slumpmässiga punkten på omkretsen projiceras sedan radiellt det erforderliga slumpmässiga avståndet med hjälp av

${\sqrt {-2\ln(s)}},\,$

använder samma s eftersom det s är oberoende av den slumpmässiga punkten på omkretsen och själv är jämnt fördelad från 0 till 1.

Genomförande

Enkel implementering i Java med medelvärde och standardavvikelse:

   
     

        
      
          
             
      
           
         
                  
                  
                    
                0
              
            
          
               
    
 privat  statisk  dubbel  reservdel  ;  privat  statisk  boolesk  hasSpare  =  false  ;  offentlig  statisk  synkroniserad  dubbel  genereraGaussian  (  dubbel  medelvärde  ,  dubbel  stdDev  )  {  if  (  hasSpare  )  {  hasSpare  =  false  ;  returnera  reserv  *  stdDev  +  medelvärde  ;  }  annat  {  dubbel  u  ,  v  ,  s  ;  gör  {  u  =  Math  .  slumpmässig  ()  *  2  -  1  ;  v  =  matematik  .  slumpmässig  ()  *  2  -  1  ;  s  =  u  *  u  +  v  *  v  ;  }  while  (  s  >=  1  ||  s  ==  );  s  =  matematik  .  sqrt  (  -2.0  *  Math  .  log  (  s  )  /  s  )  ;  reserv  =  v  *  s  ;  hasSpare  =  sant  ;  returnera  medelvärde  +  stdDev  *  u  *  s  ;  }  } 

En icke- trådssäker implementering i C++ med medelvärde och standardavvikelse:

     
      
        

      
          
             
      
           
         
                    
                    
                    
                
              
            
          
               
    
 dubbel  genereraGaussian  (  dubbelt  medelvärde  ,  dubbel  stdDev  )  {  statisk  dubbel  reservdel  ;  static  bool  hasSpare  =  false  ;  if  (  hasSpare  )  {  hasSpare  =  false  ;  returnera  reserv  *  stdDev  +  medelvärde  ;  }  annat  {  dubbel  u  ,  v  ,  s  ;  do  {  u  =  (  rand  ()  /  ((  double  )  RAND_MAX  ))  *  2,0  -  1,0  ;  v  =  (  rand  ()  /  ((  dubbel  )  RAND_MAX  ))  *  2,0  -  1,0  ;  s  =  u  *  u  +  v  *  v  ;  }  while  (  s  >=  1,0  ||  s  ==  0,0  );  s  =  sqrt  (  -2,0  *  log  (  s  )  /  s  );  reserv  =  v  *  s  ;  hasSpare  =  sant  ;  returnera  medelvärde  +  stdDev  *  u  *  s  ;  }  } 

C++11 GNU GCC libstdc++ :s implementering av std::normal_distribution använder den polära Marsaglia-metoden, som citeras här .