Newer
Older
Lecture_repo / Lectures_my / EMPP / 2016 / Lecture1 / mchrzasz.tex
@mchrzasz mchrzasz on 4 Oct 2016 34 KB Lecture 1 finished
  1. \documentclass[11 pt,xcolor={dvipsnames,svgnames,x11names,table}]{beamer}
  2.  
  3. \usepackage[english]{babel}
  4. \usepackage{polski}
  5.  
  6.  
  7. \usetheme[
  8. bullet=circle, % Other option: square
  9. bigpagenumber, % circled page number on lower right
  10. topline=true, % colored bar at the top of the frame
  11. shadow=false, % Shading for beamer blocks
  12. watermark=BG_lower, % png file for the watermark
  13. ]{Flip}
  14.  
  15. %\logo{\kern+1.em\includegraphics[height=1cm]{SHiP-3_LightCharcoal}}
  16.  
  17. \usepackage[lf]{berenis}
  18. \usepackage[LY1]{fontenc}
  19. \usepackage[utf8]{inputenc}
  20.  
  21. \usepackage{emerald}
  22. \usefonttheme{professionalfonts}
  23. \usepackage[no-math]{fontspec}
  24. \defaultfontfeatures{Mapping=tex-text} % This seems to be important for mapping glyphs properly
  25.  
  26. \setmainfont{Gillius ADF} % Beamer ignores "main font" in favor of sans font
  27. \setsansfont{Gillius ADF} % This is the font that beamer will use by default
  28. % \setmainfont{Gill Sans Light} % Prettier, but harder to read
  29.  
  30. \setbeamerfont{title}{family=\fontspec{Gillius ADF}}
  31.  
  32. \input t1augie.fd
  33.  
  34. %\newcommand{\handwriting}{\fontspec{augie}} % From Emerald City, free font
  35. %\newcommand{\handwriting}{\usefont{T1}{fau}{m}{n}} % From Emerald City, free font
  36. % \newcommand{\handwriting}{} % If you prefer no special handwriting font or don't have augie
  37.  
  38. %% Gill Sans doesn't look very nice when boldfaced
  39. %% This is a hack to use Helvetica instead
  40. %% Usage: \textbf{\forbold some stuff}
  41. %\newcommand{\forbold}{\fontspec{Arial}}
  42.  
  43. \usepackage{graphicx}
  44. \usepackage[export]{adjustbox}
  45.  
  46. \usepackage{amsmath}
  47. \usepackage{amsfonts}
  48. \usepackage{amssymb}
  49. \usepackage{bm}
  50. \usepackage{colortbl}
  51. \usepackage{mathrsfs} % For Weinberg-esque letters
  52. \usepackage{cancel} % For "SUSY-breaking" symbol
  53. \usepackage{slashed} % for slashed characters in math mode
  54. \usepackage{bbm} % for \mathbbm{1} (unit matrix)
  55. \usepackage{amsthm} % For theorem environment
  56. \usepackage{multirow} % For multi row cells in table
  57. \usepackage{arydshln} % For dashed lines in arrays and tables
  58. \usepackage{siunitx}
  59. \usepackage{xhfill}
  60. \usepackage{grffile}
  61. \usepackage{textpos}
  62. \usepackage{subfigure}
  63. \usepackage{tikz}
  64. \usepackage{hyperref}
  65. %\usepackage{hepparticles}
  66. \usepackage[italic]{hepparticles}
  67.  
  68. \usepackage{hepnicenames}
  69.  
  70. % Drawing a line
  71. \tikzstyle{lw} = [line width=20pt]
  72. \newcommand{\topline}{%
  73. \tikz[remember picture,overlay] {%
  74. \draw[crimsonred] ([yshift=-23.5pt]current page.north west)
  75. -- ([yshift=-23.5pt,xshift=\paperwidth]current page.north west);}}
  76.  
  77.  
  78.  
  79. % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
  80. \usepackage{tikzfeynman} % For Feynman diagrams
  81. \usetikzlibrary{arrows,shapes}
  82. \usetikzlibrary{trees}
  83. \usetikzlibrary{matrix,arrows} % For commutative diagram
  84. % http://www.felixl.de/commu.pdf
  85. \usetikzlibrary{positioning} % For "above of=" commands
  86. \usetikzlibrary{calc,through} % For coordinates
  87. \usetikzlibrary{decorations.pathreplacing} % For curly braces
  88. % http://www.math.ucla.edu/~getreuer/tikz.html
  89. \usepackage{pgffor} % For repeating patterns
  90.  
  91. \usetikzlibrary{decorations.pathmorphing} % For Feynman Diagrams
  92. \usetikzlibrary{decorations.markings}
  93. \tikzset{
  94. % >=stealth', %% Uncomment for more conventional arrows
  95. vector/.style={decorate, decoration={snake}, draw},
  96. provector/.style={decorate, decoration={snake,amplitude=2.5pt}, draw},
  97. antivector/.style={decorate, decoration={snake,amplitude=-2.5pt}, draw},
  98. fermion/.style={draw=gray, postaction={decorate},
  99. decoration={markings,mark=at position .55 with {\arrow[draw=gray]{>}}}},
  100. fermionbar/.style={draw=gray, postaction={decorate},
  101. decoration={markings,mark=at position .55 with {\arrow[draw=gray]{<}}}},
  102. fermionnoarrow/.style={draw=gray},
  103. gluon/.style={decorate, draw=black,
  104. decoration={coil,amplitude=4pt, segment length=5pt}},
  105. scalar/.style={dashed,draw=black, postaction={decorate},
  106. decoration={markings,mark=at position .55 with {\arrow[draw=black]{>}}}},
  107. scalarbar/.style={dashed,draw=black, postaction={decorate},
  108. decoration={markings,mark=at position .55 with {\arrow[draw=black]{<}}}},
  109. scalarnoarrow/.style={dashed,draw=black},
  110. electron/.style={draw=black, postaction={decorate},
  111. decoration={markings,mark=at position .55 with {\arrow[draw=black]{>}}}},
  112. bigvector/.style={decorate, decoration={snake,amplitude=4pt}, draw},
  113. }
  114.  
  115. % TIKZ - for block diagrams,
  116. % from http://www.texample.net/tikz/examples/control-system-principles/
  117. % \usetikzlibrary{shapes,arrows}
  118. \tikzstyle{block} = [draw, rectangle,
  119. minimum height=3em, minimum width=6em]
  120.  
  121.  
  122.  
  123.  
  124. \usetikzlibrary{backgrounds}
  125. \usetikzlibrary{mindmap,trees} % For mind map
  126. \newcommand{\degree}{\ensuremath{^\circ}}
  127. \newcommand{\E}{\mathrm{E}}
  128. \newcommand{\Var}{\mathrm{Var}}
  129. \newcommand{\Cov}{\mathrm{Cov}}
  130. \newcommand\Ts{\rule{0pt}{2.6ex}} % Top strut
  131. \newcommand\Bs{\rule[-1.2ex]{0pt}{0pt}} % Bottom strut
  132.  
  133. \graphicspath{{images/}} % Put all images in this directory. Avoids clutter.
  134.  
  135. % SOME COMMANDS THAT I FIND HANDY
  136. % \renewcommand{\tilde}{\widetilde} % dinky tildes look silly, dosn't work with fontspec
  137. \newcommand{\comment}[1]{\textcolor{comment}{\footnotesize{#1}\normalsize}} % comment mild
  138. \newcommand{\Comment}[1]{\textcolor{Comment}{\footnotesize{#1}\normalsize}} % comment bold
  139. \newcommand{\COMMENT}[1]{\textcolor{COMMENT}{\footnotesize{#1}\normalsize}} % comment crazy bold
  140. \newcommand{\Alert}[1]{\textcolor{Alert}{#1}} % louder alert
  141. \newcommand{\ALERT}[1]{\textcolor{ALERT}{#1}} % loudest alert
  142. %% "\alert" is already a beamer pre-defined
  143. \newcommand*{\Scale}[2][4]{\scalebox{#1}{$#2$}}%
  144.  
  145. \def\Put(#1,#2)#3{\leavevmode\makebox(0,0){\put(#1,#2){#3}}}
  146.  
  147. \usepackage{gmp}
  148. \usepackage[final]{feynmp-auto}
  149.  
  150. \usepackage[backend=bibtex,style=numeric-comp,firstinits=true]{biblatex}
  151. \bibliography{bib}
  152. \setbeamertemplate{bibliography item}[text]
  153.  
  154. \makeatletter\let\frametextheight\beamer@frametextheight\makeatother
  155.  
  156. % suppress frame numbering for backup slides
  157. % you always need the appendix for this!
  158. \newcommand{\backupbegin}{
  159. \newcounter{framenumberappendix}
  160. \setcounter{framenumberappendix}{\value{framenumber}}
  161. }
  162. \newcommand{\backupend}{
  163. \addtocounter{framenumberappendix}{-\value{framenumber}}
  164. \addtocounter{framenumber}{\value{framenumberappendix}}
  165. }
  166.  
  167.  
  168. \definecolor{links}{HTML}{2A1B81}
  169. %\hypersetup{colorlinks,linkcolor=,urlcolor=links}
  170.  
  171. % For shapo's formulas:
  172. \def\lsi{\raise0.3ex\hbox{$<$\kern-0.75em\raise-1.1ex\hbox{$\sim$}}}
  173. \def\gsi{\raise0.3ex\hbox{$>$\kern-0.75em\raise-1.1ex\hbox{$\sim$}}}
  174. \newcommand{\lsim}{\mathop{\lsi}}
  175. \newcommand{\gsim}{\mathop{\gsi}}
  176. \newcommand{\wt}{\widetilde}
  177. %\newcommand{\ol}{\overline}
  178. \newcommand{\Tr}{\rm{Tr}}
  179. \newcommand{\tr}{\rm{tr}}
  180. \newcommand{\eqn}[1]{&\hspace{-0.7em}#1\hspace{-0.7em}&}
  181. \newcommand{\vev}[1]{\rm{$\langle #1 \rangle$}}
  182. \newcommand{\abs}[1]{\rm{$\left| #1 \right|$}}
  183. \newcommand{\eV}{\rm{eV}}
  184. \newcommand{\keV}{\rm{keV}}
  185. \newcommand{\GeV}{\rm{GeV}}
  186. \newcommand{\im}{\rm{Im}}
  187. \newcommand{\disp}{\displaystyle}
  188. \def\be{\begin{equation}}
  189. \def\ee{\end{equation}}
  190. \def\ba{\begin{eqnarray}}
  191. \def\ea{\end{eqnarray}}
  192. \def\d{\partial}
  193. \def\l{\left(}
  194. \def\r{\right)}
  195. \def\la{\langle}
  196. \def\ra{\rangle}
  197. \def\e{{\rm e}}
  198. \def\Br{{\rm Br}}
  199. \def\fixme{{\color{red} FIXME!}}
  200. \def\mc{{\color{Magenta}{MC}}}
  201. \def\pdf{{\rm p.d.f.}}
  202.  
  203. \author{ {\fontspec{Trebuchet MS}Marcin Chrz\k{a}szcz} (Universit\"{a}t Z\"{u}rich)}
  204. \institute{UZH}
  205. \title[Introduction to \\Monte Carlo methods]{Introduction to \\Monte Carlo methods}
  206. \date{\fixme}
  207.  
  208.  
  209. \begin{document}
  210. \tikzstyle{every picture}+=[remember picture]
  211.  
  212. {
  213. \setbeamertemplate{sidebar right}{\llap{\includegraphics[width=\paperwidth,height=\paperheight]{bubble2}}}
  214. \begin{frame}[c]%{\phantom{title page}}
  215. \begin{center}
  216. \begin{center}
  217. \begin{columns}
  218. \begin{column}{0.9\textwidth}
  219. \flushright\fontspec{Trebuchet MS}\bfseries \Huge {Introduction to \\Monte Carlo methods}
  220. \end{column}
  221. \begin{column}{0.2\textwidth}
  222. %\includegraphics[width=\textwidth]{SHiP-2}
  223. \end{column}
  224. \end{columns}
  225. \end{center}
  226. \quad
  227. \vspace{3em}
  228. \begin{columns}
  229. \begin{column}{0.44\textwidth}
  230. \flushright \vspace{-1.8em} {\fontspec{Trebuchet MS} \Large Marcin Chrząszcz\\\vspace{-0.1em}\small \href{mailto:mchrzasz@cern.ch}{mchrzasz@cern.ch}}
  231.  
  232. \end{column}
  233. \begin{column}{0.53\textwidth}
  234. \includegraphics[height=1.3cm]{uzh-transp}
  235. \end{column}
  236. \end{columns}
  237.  
  238. \vspace{1em}
  239. % \footnotesize\textcolor{gray}{With N. Serra, B. Storaci\\Thanks to the theory support from M. Shaposhnikov, D. Gorbunov}\normalsize\\
  240. \vspace{0.5em}
  241. \textcolor{normal text.fg!50!Comment}{Experimental Methods in Particle Physics, \\ 5 October, 2016}
  242. \end{center}
  243. \end{frame}
  244. }
  245.  
  246.  
  247. \begin{frame}\frametitle{Literature}
  248.  
  249. \begin{enumerate}
  250. \item J. M. Hammersley, D. C. Hamdscomb, ``Monte Carlo Methods'', London: Methuen \& Co. Ltd., New York: J. Wiley \& Sons Inc., 1964
  251. \item I. M. Sobol, ``The Monte Carlo Method'', Mir Publishers, Moscow, 1975.
  252. \item M. H. Kalos, P. A. Whitlock, ,,Monte Carlo Methods”, J. Wiley \& Sons Inc., New York, 1986
  253. \item G. S. Fishman, ,,Monte Carlo: Concepts, Algorithms and Applications”, Springer, 1996.
  254. \item R. Y. Rubinstein, D. P. Kroese, ,,Simulation and the Monte Carlo Method”, Second Edition, J. Wiley \& Sons Inc., 2008.
  255. \item R. Korn, E. Korn, G. Kroisandt, ,,Monte Carlo methods and models in finance and insurance”, CRC Press,
  256. Taylor \& Francis Group, 2010.
  257. \item S. Jadach, ,,Practical Guide to Monte Carlo”, \href{http://arxiv.org/abs/physics/9906056}{arXiv:physics/9906056}, \href{http://cern.ch/jadach/MCguide/}{http://cern.ch/jadach/MCguide/}.
  258. \end{enumerate}
  259.  
  260.  
  261. \end{frame}
  262.  
  263.  
  264.  
  265.  
  266. \begin{frame}\frametitle{Course Plan}
  267.  
  268. We will have 6 hours of Monte Carlo (MC) lectures. The lectures will be devoted:\\
  269.  
  270. \hspace{1.5cm}
  271. \begin{itemize}
  272. \item 1 h: Mathematical introduction to MC methods.
  273. \item 1 h: MC integration methods.
  274. \item 2 h: Random numbers generators.
  275. \item 2 h: Markov Chain MC.
  276. \item 2 h: Tutorial and examples.
  277. \end{itemize}
  278. \hspace{1cm} \\
  279. The hands-on tutorial will consist of program templates in which we will implement couple of algorithms that were explained in the lecture. \\
  280. $\Rrightarrow$ All examples shown in this course are available in the github repository:\\
  281. \url{https://github.com/mchrzasz/EMPP_MC}\\
  282. \color{RubineRed}{There will be an indication (in this color) on the adequate slide for each of the macro.}
  283.  
  284.  
  285. \end{frame}
  286.  
  287.  
  288. \begin{frame}\frametitle{Definitions}
  289. \begin{footnotesize}
  290. $\Rrightarrow$ Basic definition:\\
  291. \begin{exampleblock}{}
  292. Monte Carlo method is any technique that uses {\it{random numbers}} to solve a given mathematical problem.
  293. \end{exampleblock}
  294. % \vspace{0.5cm}
  295.  
  296. $\rightarrowtail$ Random number: For the purpose of this course we need to assume that we know what it is, although the formal definition is highly non-trivial.\\
  297. \vspace{0.05cm}
  298. $\Rrightarrow$ My favourite definition (Halton 1970): \begin{scriptsize}more complicated, but more accurate.\end{scriptsize}
  299.  
  300. \begin{exampleblock}{}
  301. ''Representing the solution of a problem as a parameter of a hypothetical population, and using a random sequence of numbers to construct a sample of the population, from which statistical estimates of the parameter can be obtained.''
  302. \end{exampleblock}
  303. To put this definition in mathematical language:\\
  304. Let $F$ be a solution of a given mathematical problem. The estimate of the result $\hat{F}$:\\
  305. \begin{equation*}
  306. \hat{F}=f( \lbrace r_1, r_2, r_3,...,r_n \rbrace; ...),
  307. \end{equation*}
  308. where $\lbrace r_1, r_2, r_3,...,r_n \rbrace$ are random numbers.
  309. \begin{center}
  310. \color{red}{The problem we are solving doesn't need to be stochastic!}
  311. \end{center}
  312. \begin{scriptsize}
  313. $\twoheadrightarrow$ One could wonder why are we trying to add all the stochastic properties to a deterministic problem. Those are the properties that allow to use all well known statistic theorems.
  314. \end{scriptsize}
  315. \end{footnotesize}
  316. \end{frame}
  317.  
  318. \begin{frame}\frametitle{History of MC methods}
  319. \begin{footnotesize}
  320. \begin{itemize}
  321. \item {\color{PineGreen}{G. Compte de Buffon (1777)}} - First documented usage of random numbers for integral computation (Buffon thrown niddle on the table with parrarel line; we will do a modern version of this exercise).
  322. \item {\color{PineGreen}{Marquis de Laplace (1886)}} - Used the Buffon niddle to determine the value of $\pi$ number.
  323. \item {\color{PineGreen}{Lord Kelvin (1901)}} - Thanks to drawing randomly numbered cards he managed he managed to calculate some integrals in kinematic gas theorem.
  324. \item {\color{PineGreen}{W. S. Gosse (better knows as Student) (1908)}} - Used similar way as Lord Kelvin to get random numbers to prove \textit{t}-Student distribution.
  325. \item {\color{PineGreen}{Enrico Fermi (1930) }} - First mechanical device (\texttt{FERMIAC}) for random number generations. Solved neutron transport equations in the nuclear plants.
  326. \item {\color{PineGreen}{S. Ulam, R. Feynman, J. von Neumann et. al.}} - First massive usage of random numbers. Most applications were in Manhattan project to calculate neutron scattering and absorption. \\
  327. In {\color{NavyBlue}{Los Alamos}} the name {\color{Mahogany}{Monte Carlo}} was created as kryptonim of this kind of calculations.
  328. \end{itemize}
  329.  
  330.  
  331. \end{footnotesize}
  332. \end{frame}
  333.  
  334.  
  335.  
  336. \begin{frame}\frametitle{Euler number determination, $\rm \color{RubineRed}{Lecture1/Euler\_number}$ }
  337. \begin{footnotesize}
  338. $\Rrightarrow$ As mentioned before \mc~methods can be used to solve problems that \textbf{do not} have stochastic nature! All the integrals calculated in Los Alamos during the Manhattan project are nowadays solvable without any \mc~methods.\\
  339. $\rightarrowtail$ Let's give a trivial example of solving a non stochastic problem: calculating Euler number $e$. We know that $e=2.7182818...$.
  340. $\Rrightarrow$ To calculate the $\hat{e}$ we will use the following algorithm:
  341. \begin{itemize}
  342. \item We generate a random number in range $(0,1)$ (in stat. $\mathcal{U}(0,1)$) until the number we generate is smaller then the previous one, aka we get the following sequence:
  343. \begin{align*}
  344. x_1<x_2<...<x_{n-1}>x_{n}
  345. \end{align*}
  346. \item We store the number $n$. We repeat this experiment $N$ times and calculate the arithmetic average of $n$. The obtained value is an statistical estimator of $e$:
  347. \begin{align*}
  348. \hat{e}= \dfrac{1}{N}\sum_{i=1}^N n_i \xrightarrow{N\to \infty} e .
  349. \end{align*}
  350. \end{itemize}
  351. $\Rrightarrow$ Numerical example:
  352. \begin{tabular}{r c c c }
  353. $N$ & $\hat{e}$ & $\hat{e} - e$ & \multirow{5}{*}{Is this $\sim\sqrt{N}$?} \\
  354. 100 & $2.760000$ & $0.041718$ \\
  355. 10000 & $2.725000$ & $0.006718$ \\
  356. 1000000 & $2.718891$ & $0.000609$ \\
  357. 100000000 & $2.718328$ & $0.000046$\\
  358. \end{tabular}
  359.  
  360.  
  361. \end{footnotesize}
  362. \end{frame}
  363.  
  364. \begin{frame}\frametitle{Let's test the $\sqrt{N}$, $\rm \color{RubineRed}{Lecture1/Euler\_number}$ }
  365.  
  366. \only<1>
  367. {
  368. $\Rrightarrow$ In the last example we measured the Euler number using different number of pseudo-experiments.\\
  369. $\rightarrowtail$ We compared the obtained value to the true and observed roughly a $\sqrt{N}$ dependence on the difference between the true value and the obtained one.\\
  370. }
  371. $\rightarrowtail$ Could we test this? YES! Lets put our experimentalist hat on!\\
  372. $\rightarrowtail$ From the begging of studies they tooth us to get the error you need to repeat the measurements.
  373. \only<1>
  374. {
  375. \begin{exampleblock}{The algorithm:}
  376. Previous time we measured Euler number using $N$ events, where $N \in (100, 1000, 10000, 100000)$. Now lets repeat this measurement $n_N$ times (of course each time we use new generated numbers). From the distribution of $\hat{e} -e$ we could say something about the uncertainty of our estimator for given $N$.
  377. \end{exampleblock}
  378. }
  379.  
  380. \begin{center}
  381. \only<2>
  382. {
  383. \includegraphics[angle=-90,width=0.8\textwidth]{images/result_error.pdf}
  384. }
  385. \only<3>
  386. {
  387. \includegraphics[angle=-90,width=0.8\textwidth]{images/result_error_dep.pdf}
  388. }
  389. \end{center}
  390.  
  391.  
  392.  
  393.  
  394.  
  395.  
  396.  
  397.  
  398.  
  399.  
  400.  
  401.  
  402. \end{frame}
  403.  
  404. \begin{frame}\frametitle{Monte Carlo and integration}
  405. \begin{footnotesize}
  406. $\hookrightarrow$ {\color{BrickRed}{\textbf{All MC calculations are equivalent to preforming an integration.}}}\\
  407. $\rightrightarrows$ Assumptions: $r_i$ random numbers from $\mathcal{U}(0,1)$. The MC result:
  408. \begin{align*}
  409. F=F(r_1,r_2,...r_n)
  410. \end{align*}
  411. is unbias estimator of an integral:
  412. \begin{align*}
  413. I=\int_0^1...\int_0^1 F(x_1,x_2,...,x_n)dx_1,dx_2...,dx_n
  414. \end{align*}
  415. aka the expected value of the $I$ integral is:
  416. \begin{align*}
  417. E(F)=I.
  418. \end{align*}
  419. \begin{exampleblock}{}
  420. $\Rrightarrow$ This mathematical identity is the most useful property of the MC methods. It is a link between mathematical analysis and statistic world. Now we can use the best of the both world!
  421. \end{exampleblock}
  422. If we want to calculate the integral in different range then $(0,1)$ we just scale the the previous result:
  423. \begin{align*}
  424. \dfrac{1}{N}\sum_{i=1}^N f(x_i) \xrightarrow{N\to \infty} E(f)=\dfrac{1}{b-a}\int_a^b f(x)dx
  425. \end{align*}
  426.  
  427.  
  428. \end{footnotesize}
  429. \end{frame}
  430.  
  431.  
  432. \begin{frame}\frametitle{Uncertainty from Monte Carlo methods}
  433. \begin{footnotesize}
  434. $\Rrightarrow$ In practice we do not have $N\to \infty$ so we will never know the exact result of an integral :(\\
  435. $\longmapsto$ Let's use the {\color{BrickRed}{statistical}} world and estimate the uncertainty of an integral in this case :)\\
  436. $\rightarrowtail$ A variance of a MC integral:
  437. \begin{align*}
  438. V(\hat{I}) = \dfrac{1}{n} \Big\lbrace E(f^2) - E^2(f) \Big\rbrace = \dfrac{1}{n} \Big\lbrace \dfrac{1}{b-a} \int_a^b f^2(x)dx - I^2 \Big\rbrace
  439. \end{align*}
  440. \begin{alertblock}{}
  441. $\looparrowright$ To calculate $V(\hat{I})$ one needs to know the value of $I$!
  442. \end{alertblock}
  443. $\Rrightarrow$ In practice $V(\hat{I})$ is calculated via estimator:
  444. \begin{columns}
  445. \column{2in}
  446. \begin{align*}
  447. \hat{V}(\hat{I})=\dfrac{1}{n}\hat{V}(f),
  448. \end{align*}
  449. \column{3in}
  450. \begin{align*}
  451. \hat{V}(f) = \dfrac{1}{n-1}\sum_{i=1}^n \Big[ f(x_i)-\dfrac{1}{n} \sum_{i=1}^nf(x_i)\Big]^2.
  452. \end{align*}
  453. \end{columns}
  454.  
  455.  
  456. $\Rrightarrow$ MC estimator of standard deviation: $\hat{\sigma}=\sqrt{\hat{V}(\hat{I})}$
  457.  
  458.  
  459. \end{footnotesize}
  460. \end{frame}
  461.  
  462.  
  463. \begin{frame}\frametitle{Buffon needle - $\pi$ number calculus}
  464. \begin{footnotesize}
  465.  
  466. $\Rrightarrow$ Buffon needle (Buffon 1777, Laplace 1886):
  467. We are throwing a needle (of length $l$) on to a surface covered with parallel lines width distance $L$. If a thrown needle touches a line we count a hit, else miss. Knowing the number of hits and misses one can calculate the $\pi$ number.
  468. \vspace{0.3cm}
  469. \begin{columns}
  470. \column{0.1in}
  471. {~}
  472. \column{2in}
  473. {\color{ForestGreen}{Experiment:}}
  474. \column{2.8in}
  475. {\color{Cerulean}{Theory:}}
  476.  
  477. \end{columns}
  478.  
  479.  
  480. \begin{columns}
  481. \column{0.1in}
  482. {~}
  483. \column{2in}
  484.  
  485. \includegraphics[width=0.9\textwidth]{images/buffon.png}\\
  486. $n$ - number of hits\\
  487. $N$ number of hits and misses,\\
  488. aka number of tries.
  489.  
  490. \column{2.8in}
  491. $\Rightarrow$ x - angle between needle and horizontal line, $x \in \mathcal{U}(0,\pi)$.
  492. $\Rightarrow$ the probability density function (\pdf) for x:
  493. \begin{align*}
  494. \rho(x)=\dfrac{1}{\pi}
  495. \end{align*}
  496. $\Rightarrow$ $p(x)$ probability to hit a line for a given x value:
  497. \begin{align*}
  498. p(x)=\dfrac{l}{L}\vert \cos x \vert
  499. \end{align*}
  500. $\Rightarrow$ Total hit probability:
  501. \begin{align*}
  502. P = E[p(x)]=\int_0^{\pi}p(x)\rho(x)dx=\dfrac{2l}{\pi L}
  503. \end{align*}
  504.  
  505. \end{columns}
  506. Now one can calculate $\hat{P}$ from MC : $\hat{P}=\dfrac{n}{N} \xrightarrow{N\to \infty} P= \dfrac{2l}{\pi L} \Rightarrow \hat{\pi}=\dfrac{2Nl}{nL}$
  507.  
  508.  
  509. \end{footnotesize}
  510. \end{frame}
  511.  
  512.  
  513.  
  514.  
  515.  
  516.  
  517.  
  518.  
  519.  
  520.  
  521.  
  522. \begin{frame}\frametitle{Buffon needle - Simplest Carlo method}
  523. \begin{footnotesize}
  524. {\color{MidnightBlue}{Monte Carlo type ''heads or tails''}}\\
  525. Let's use the summery of $p(x)$ function nad take $0<x<\frac{\pi}{2}$.\\
  526. $\Rightarrow$ Algorithm:\\
  527. \begin{columns}
  528. \column{0.1in}
  529. {~}
  530. \column{3.2in}
  531.  
  532.  
  533.  
  534. Generate 2 dim. distribution:
  535. \begin{align*}
  536. (x,y): \mathcal{U}(0,\dfrac{\pi}{2})\times \mathcal{U}(0,1) {\rm{~and~}}
  537. \end{align*}
  538. \begin{align*}
  539. y
  540. \begin{cases}
  541. \leq p(x): & \text{hit},\\
  542. > p(x): & \text{miss}.
  543. \end{cases}
  544. \end{align*}
  545.  
  546. \column{2.5in}
  547. \includegraphics[width=0.75\textwidth]{images/result.png}
  548.  
  549.  
  550.  
  551. \end{columns}
  552. Let's define weight function: $w(x,y)=\Theta(p(x)-y)$, \\
  553. where $\Theta(x)$ is the step function.\\
  554. $\rightarrowtail$ \pdf : $\varrho(x,y)=\rho(x)g(y)=\frac{2}{\pi} \cdot 1$\\
  555. $\Rightarrow$ Integrated probability:
  556. \begin{align*}
  557. P=E(w)=\int w(x,y) \varrho(x,y)dx dy = \dfrac{2l}{\pi L} \xleftarrow{N\to \infty}\hat{P}=\frac{1}{N} \sum_{i=1}^N w(x_i,y_i)= \dfrac{n}{N}
  558. \end{align*}
  559. Standard deviation for $\hat{P}$: $\hat{\sigma}=\dfrac{1}{\sqrt{N-1}}\sqrt{\dfrac{n}{N}\Big(1-\dfrac{n}{N}\Big)} $
  560.  
  561.  
  562.  
  563. \end{footnotesize}
  564. \end{frame}
  565.  
  566.  
  567. \begin{frame}\frametitle{Buffon needle, $\rm \color{RubineRed}{Lecture1/Heads\_tails}$ }
  568. \begin{small}
  569.  
  570.  
  571. $\Rrightarrow$ Lets make this toy experiment and calculate the $\pi$ number.\\
  572. $\hookrightarrow$ We can simulate the central position $(y)$ of an needle between $(-L, L)$ from $\mathcal{U}(-L, L)$.
  573.  
  574. \begin{exampleblock}{Symmetry:}
  575. Please note the symmetry of the problem, if the position of the needle would be $>L$ then we can shift the needle by any number of $L$'s.
  576. \end{exampleblock}
  577.  
  578. $\hookrightarrow$ New we simulate the angle $(\phi)$ with a flat distribution from $(0,\pi)$.
  579. $\hookrightarrow$ The maximum and minimum $y$ position of the needle are:
  580. \begin{align*}
  581. y_{\max}=y+\vert \cos \phi \vert l\\
  582. y_{\min}=y-\vert \cos \phi \vert l
  583. \end{align*}
  584. $\hookrightarrow$ Now we check if the needle touches any of the lines: $y=L$, $y=0$ or $y=-L$. If yes we count the events.
  585. \end{small}
  586. \begin{center}
  587.  
  588.  
  589. \begin{footnotesize}
  590.  
  591.  
  592. \begin{tabular}{r c c c }
  593. $N$ & $\hat{\pi}$ & $\hat{\pi} - \pi$ & $\sigma(\hat{\pi})$ \\
  594. 10000 & $3.12317$ & $-0.01842$ & $0.03047$\\
  595. 100000 & $3.14707$ & $0.00547$ & $0.00979$\\
  596. 1000000 & $3.13682$ & $-0.00477$ & $0.00307$\\
  597. 10000000 & $3.14096$ & $-0.00063$ & $0.00097$\\
  598.  
  599.  
  600. \end{tabular}
  601.  
  602. \end{footnotesize}
  603. \end{center}
  604. \end{frame}
  605. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  606.  
  607. \begin{frame}\frametitle{Central Limit Theorem, $\rm \color{RubineRed}{Lecture1/CLT}$ }
  608. \begin{footnotesize}
  609. \begin{exampleblock}{}
  610. Large independent random numbers assembly has always Gaussian distribution no matter from what distribution they were generated from as far as they have finite variances and expected values and the assembly is sufficiently large.
  611. \end{exampleblock}
  612. \includegraphics[width=0.9\textwidth]{images/dupa.png}
  613.  
  614.  
  615. \end{footnotesize}
  616. \end{frame}
  617.  
  618.  
  619. \begin{frame}\frametitle{Crude Monte Carlo method of integration}
  620. \begin{footnotesize}
  621. $\Rrightarrow$ {\color{MidnightBlue}{Crude Monte Carlo method of integration is based on Central Limit Theorem (CLT): }}\\
  622. \begin{align*}
  623. \dfrac{1}{N} \sum_{i=1}^N f(x_i) \xrightarrow{N\to \infty} \dfrac{1}{b-a}\int_a^b f(x)dx =E(f)
  624. \end{align*}
  625. $\Rrightarrow$ The standard deviation can be calculated:
  626. \begin{align*}
  627. \sigma = \dfrac{1}{\sqrt{N}} \sqrt{\Big[ E(f^2) -E^2(f)\Big] }
  628. \end{align*}
  629.  
  630. $\Rrightarrow$ From LNT we have:
  631. \begin{align*}
  632. P= \int w(x) \rho(x) dx = \int_0^{\pi/2} (\frac{l}{L} \cos x ) \frac{2}{\pi} dx= \dfrac{2l}{\pi L} \xrightarrow{N\to \infty} \frac{1}{N}\sum_{i=1}^N w(x_i)
  633. \end{align*}
  634. $\Rrightarrow$ Important comparison between ''Hit and mishit'' and Crude \mc~methods. One can analytically calculate:
  635.  
  636. \begin{align*}
  637. \hat{\sigma}^{{\rm{Crude}}} < \hat{\sigma}^{{\rm{Hit~and~mishit}}}
  638. \end{align*}
  639.  
  640.  
  641. $\Rrightarrow$ Crude \mc~is \textbf{always} better then ''Hit and mishit'' method. We will prove this on an example (can be proven analytically as well).
  642.  
  643.  
  644. \end{footnotesize}
  645. \end{frame}
  646.  
  647. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  648. \begin{frame}\frametitle{Crude MC vs Hit and mishit”, $\rm \color{RubineRed}{Lecture1/Crude\_vs\_HT}$}
  649.  
  650. $\Rrightarrow$ We can repeat a toy MC studies as we did in the Euler needle case.\\
  651. $\hookrightarrow$ In this example we want to calculate $\int_0^{\pi/2} \cos x dx$
  652.  
  653. \only<1>{
  654. \begin{center}
  655. \includegraphics[angle=-90,width=0.8\textwidth]{images/results_0.pdf}
  656. \end{center}
  657. }
  658. \only<2>{
  659. \begin{center}
  660. \includegraphics[angle=-90,width=0.8\textwidth]{images/results_1.pdf}
  661. \end{center}
  662. }
  663. \only<3>{
  664. \begin{center}
  665.  
  666. \includegraphics[angle=-90,width=0.4\textwidth]{images/results_fit_0.pdf}
  667. \includegraphics[angle=-90,width=0.4\textwidth]{images/results_fit_1.pdf}\\
  668. \end{center}
  669. $\Rrightarrow$ One clearly sees that both methods follow $1/\sqrt{N}$ dependence and that the Crude MC is always better then the ''Hit and mishit''.\\
  670. $\Rrightarrow$ Please note that for the ''Hit and mishit'' we are suing 2 times more random numbers than for the Crude method so in terms of timing the Crude MC is also much faster.
  671.  
  672. }
  673.  
  674.  
  675.  
  676. \end{frame}
  677.  
  678. \begin{frame}\frametitle{Classical methods of variance reduction}
  679. \begin{footnotesize}
  680.  
  681. $\Rrightarrow$ In Monte Carlo methods the statistical uncertainty is defined as:
  682. \begin{align*}
  683. \sigma = \dfrac{1}{\sqrt{N}}\sqrt{V(f)}
  684. \end{align*}
  685. $\Rrightarrow$ Obvious conclusion:
  686. \begin{itemize}
  687. \item To reduce the uncertainty one needs to increase $N$.\\
  688. $\rightrightarrows$ Slow convergence. In order to reduce the error by factor of 10 one needs to simulate factor of 100 more points!
  689. \end{itemize}
  690. $\Rrightarrow$ How ever the other handle ($V(f)$) can be changed! $\longrightarrow$ Lot's of theoretical effort goes into reducing this factor.\\
  691. $\Rrightarrow$ We will discuss {\color{Mahogany}{four}} classical methods of variance reduction:
  692. \begin{enumerate}
  693. \item Stratified sampling.
  694. \item Importance sampling.
  695. \item Control variates.
  696. \item Antithetic variates.
  697. \end{enumerate}
  698.  
  699.  
  700.  
  701.  
  702.  
  703. \end{footnotesize}
  704. \end{frame}
  705.  
  706.  
  707.  
  708.  
  709. \begin{frame}\frametitle{Stratified sampling}
  710. \begin{footnotesize}
  711. $\Rrightarrow$ The most intuitive method of variance reduction. The idea behind it is to divide the function in different ranges and to use the Riemann integral property:
  712. \begin{align*}
  713. I = \int_0^1 f(u) du = \int_0^a f(u)du + \int_a^1 f(u) du,~ 0<a<1.
  714. \end{align*}
  715.  
  716. $\Rrightarrow$ The reason for this method is that in smaller ranges the integration function is more flat. And it's trivial to see that the more flatter you get the smaller uncertainty.
  717. $\rightrightarrows$ A constant function would have zero uncertainty!
  718.  
  719. \begin{exampleblock}{General schematic:}
  720. Let's take our integration domain and divide it in smaller domains. In the $j^{th}$ domain with the volume $w_j$ we simulate $n_j$ points from uniform distribution. We sum the function values in each of the simulated points for each of the domain. Finally we sum them with weights proportional to $w_i$ and anti-proportional to $n_i$.
  721. \end{exampleblock}
  722.  
  723. \end{footnotesize}
  724. \end{frame}
  725.  
  726. \begin{frame}\frametitle{Stratified sampling - mathematical details}
  727. \begin{footnotesize}
  728. Let's define our integrals and domains:
  729. \begin{align*}
  730. I=\int_{\Omega} f(x) dx,{~}{~} \Omega=\bigcup_{i=1}^k w_i
  731. \end{align*}
  732. The integral over $j^{th}$ domain:
  733. \begin{align*}
  734. I_j=\int_{w_j} f(x) dx,{~}{~} \Rightarrow I = \sum_{j=1}^k I_i
  735. \end{align*}
  736. $\rightrightarrows$ $p_j$ uniform distribution in the $w_j$ domain: $dp_j=\frac{dx}{w_j}$.\\
  737. $\rightrightarrows$ The integral is calculated based on crude \mc~method. The estimator is equal:
  738. \begin{align*}
  739. \hat{I}_j = \frac{w_j}{n_j} \sum_{i=1}^{n_j}f(x_j^i)
  740. \end{align*}
  741. Now the total integral is just a sum:
  742. \begin{align*}
  743. \hat{I} = \sum_{j=1}^k \hat{I}_j = \sum_{j=1}^k \frac{w_j}{n_j} \sum_{i=1}^{n_j} f(x_j^{(i)}),
  744. \end{align*}
  745. \begin{columns}
  746. \column{0.2in}
  747. \column{2in}
  748. Variance:
  749. $V(\hat{I})=\sum_{j=1}^k \dfrac{w_j^2}{n_j}V_j(f)$,
  750. \column{3in}
  751. and it's estimator:
  752. $\hat{V} (\hat{I}) = \sum_{j=1}^k \dfrac{w_j^2}{n_j} \hat{V}_j(f)$
  753. \end{columns}
  754.  
  755. \end{footnotesize}
  756. \end{frame}
  757.  
  758. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  759.  
  760. \begin{frame}\frametitle{Importance sampling}
  761. \begin{footnotesize}
  762. $\Rrightarrow$ If the function is changing rapidly in its domain one needs to use a more elegant method: make the function more stable.\\
  763. $\rightrightarrows$ The solution is from first course of mathematical analysis: change the integration variable :)
  764. \begin{align*}
  765. f(x)dx \longrightarrow \frac{f(x)}{g(x)}dG(x),~{\rm{where}}~g(x)=\frac{dG(x)}{dx}
  766. \end{align*}
  767. \begin{exampleblock}{Schematic:}
  768. \begin{itemize}
  769. \item Generate the distribution from $G(x)$ instead of $\mathcal{U}$.
  770. \item For each generate point calculate the weight: $w(x)=\frac{f(x)}{g(x)}$.
  771. \item We calculate the expected value $\hat{E}(w)$ and its variance $\hat{V}_G(w)$ for the whole sample.
  772. \end{itemize}
  773. \end{exampleblock}
  774.  
  775.  
  776.  
  777. \begin{itemize}
  778.  
  779. \item If $g(x)$ is choose correctly the resulting variance can be much smaller.
  780. \item There are some mathematical requirements:
  781. \begin{itemize}
  782. \item $g(x)$ needs to be non-negative and analytically integrable on its domain.
  783. \item $G(x)$ invertible or there should be a direct generator of $g$ distribution.
  784. \end{itemize}
  785. \end{itemize}
  786.  
  787.  
  788. \end{footnotesize}
  789. \end{frame}
  790.  
  791.  
  792. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  793.  
  794. \begin{frame}\frametitle{Importance sampling - Example}
  795. \begin{footnotesize}
  796.  
  797. $\Rrightarrow$ Let's take our good old $\pi$ determination example.
  798. \vspace{0.3cm}
  799. \begin{columns}
  800. \column{0.1in}
  801. \column{3in}
  802. $\Rrightarrow$ Let's take here for simplicity: $L=l$.
  803. \begin{itemize}
  804. \item Let's take a trivial linear weight function: $g(x)=\frac{4}{\pi}(1-\frac{2}{\pi}x)$
  805. \item It's invertible analytically: $G(x)=\frac{4}{\pi}x(1-\frac{x}{\pi})$
  806. \item The weight function:
  807. \begin{align*}
  808. w(x)=\frac{p(x)}{g(x)}=\frac{\pi}{4}\frac{\cos x}{1-2x/ \pi}
  809. \end{align*}
  810. \end{itemize}
  811.  
  812. \column{2in}
  813. \includegraphics[width=0.95\textwidth]{images/result_weight.png}
  814. \end{columns}
  815. \begin{itemize}
  816. \item Now the new standard deviation is smaller:
  817. \end{itemize}
  818. \begin{align*}
  819. \sigma_{\pi}^{{\rm{IS}}} \simeq \frac{0.41}{\sqrt{N}} < \sigma_{\pi}\simeq \frac{1.52}{\sqrt{N}}
  820. \end{align*}
  821. \begin{itemize}
  822. \item Importance sampling has advantages:
  823. \begin{itemize}
  824. \item Big improvements of variance reduction.
  825. \item The only method that can cope with singularities.
  826. \end{itemize}
  827. \end{itemize}
  828.  
  829.  
  830. \end{footnotesize}
  831. \end{frame}
  832.  
  833.  
  834. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  835.  
  836.  
  837. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  838.  
  839. \begin{frame}\frametitle{Wrap up}
  840. \begin{footnotesize}
  841. $\Rrightarrow$ To sum up:
  842. \begin{itemize}
  843. \item We discussed basic mathematical properties of \mc~methods.
  844. \item We shown that besides the stochastic nature of \mc~ they can be used to determine totally non stochastic quantities.
  845. \item We demonstrated there is a perfect isomorphism between \mc~method and integration.
  846. \item We learned how co calculate integrals and estimate the uncertainties.
  847. \item Finally we discussed several classical methods of variance reduction.
  848. \end{itemize}
  849.  
  850.  
  851. \end{footnotesize}
  852. \end{frame}
  853.  
  854.  
  855.  
  856.  
  857.  
  858.  
  859.  
  860.  
  861. \backupbegin
  862.  
  863. \begin{frame}\frametitle{Backup}
  864.  
  865.  
  866. \end{frame}
  867.  
  868. \begin{frame}\frametitle{Control variates}
  869. \begin{footnotesize}
  870. $\Rrightarrow$ Control variates uses an other nice property of Riemann integral:
  871. \begin{align*}
  872. \int f(x) dx = \int [f(x)-g(x)]dx+ \int g(x)dx
  873. \end{align*}
  874. \begin{itemize}
  875. \item $g(x)$ needs to be analytically integrable.
  876. \item The uncertainty comes only from the integral: $\int [f(x)-g(x)]dx$.
  877. \item Obviously: $V(f\to g) \xrightarrow{f\to g} 0$
  878. \end{itemize}
  879. $\Rrightarrow$ Advantages:\\
  880. \begin{itemize}
  881. \item Quite stable, immune to the singularities.
  882. \item $g(x)$ doesn't need to be invertible analytically.
  883. \end{itemize}
  884. $\Rrightarrow$ Disadvantage:\\
  885. \begin{itemize}
  886. \item Useful only if you know $\int g(x)dx$
  887. \end{itemize}
  888.  
  889. \end{footnotesize}
  890. \end{frame}
  891.  
  892.  
  893. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  894.  
  895. \begin{frame}\frametitle{Antithetic variates}
  896. \begin{footnotesize}
  897. $\Rrightarrow$ In \mc~methods usually one uses the independent random variables. The Antithetic variates method on purpose uses a set of correlated variables (negative correlation is the important property):
  898. \begin{itemize}
  899. \item Let $f$ and $f\prime$ will be functions of x on the same domain.
  900. \item The variance: $V(f+f\prime)=V(f)+V(f\prime)+2 Cov(f,f\prime)$.
  901. \item If $Cov(f,f\prime)<0$ then you can reduce the variance.
  902. \end{itemize}
  903. $\Rrightarrow$ Advantages:
  904. \begin{itemize}
  905. \item If you can pick up $f$ and $f\prime$ so that they have negative correlation one can significantly reduce the variance!
  906. \end{itemize}
  907.  
  908. $\Rrightarrow$ Disadvantages:
  909. \begin{itemize}
  910. \item There are no general methods to produce such a negative correlations.
  911. \item Hard to generalize this for multidimensional case.
  912. \item You can't generate events from $f(x)$ with this method.
  913. \end{itemize}
  914.  
  915.  
  916.  
  917. \end{footnotesize}
  918. \end{frame}
  919.  
  920. \backupend
  921.  
  922. \end{document}