<p><div><div><div><div data-show-faces="false" data-height="21" data-width="100" data-layout="button_count" data-send="true" data-href="https://www.r-bloggers.com/uniform-on-the-sphere-or-not/"></div></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="https://xianblog.wordpress.com/2018/03/08/uniform-on-the-sphere-or-not/"> R – Xi'an's Og</a></strong>, and kindly contributed to <a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers)</a> </div>
<p><strong><a rel="nofollow" href="https://xianblog.files.wordpress.com/2018/03/sinsim.jpg" target="_blank"><img data-orig-file="https://xianblog.files.wordpress.com/2018/03/sinsim.jpg?w=298&amp;h=298&amp;fit=298%2C298&amp;resize=298%2C298" data-orig-size="480,480" data-medium-file="https://xianblog.files.wordpress.com/2018/03/sinsim.jpg?w=298&amp;h=298&amp;fit=298%2C298&amp;resize=298%2C298?w=300" data-attachment-id="35584" data-image-meta="{" data-image-description="" src="https://xianblog.files.wordpress.com/2018/03/sinsim.jpg?w=298&h=298&fit=298%2C298&resize=298%2C298" data-permalink="https://xianblog.wordpress.com/2018/03/08/uniform-on-the-sphere-or-not/sinsim/" data-large-file="https://xianblog.files.wordpress.com/2018/03/sinsim.jpg?w=298&amp;h=298&amp;fit=298%2C298&amp;resize=298%2C298?w=450" alt="" data-image-title="sinsim" sizes="(max-width: 298px) 100vw, 298px" width="298" data-comments-opened="1" class="wp-image-35584" data-recalc-dims="1" height="298"></a>W</strong>hile looking at <a rel="nofollow" href="https://stats.stackexchange.com/q/331253/7224" target="_blank">X validated questions</a>, I came upon <a rel="nofollow" href="https://stats.stackexchange.com/q/331253/7224" target="_blank">this comment</a> that simulating a uniform distribution on a d-dimensional unit sphere does not proceed from generating angles at random on (0,2π) and <a rel="nofollow" href="https://stats.stackexchange.com/q/77440/7224" target="_blank">computing spherical coordinates</a>… Which I must confess would have been my initial suggestion! This is obvious, nonetheless, when computing <a rel="nofollow" href="http://corysimon.github.io/articles/uniformdistn-on-sphere/" target="_blank">the Jacobian</a> of the spherical coordinate transform, which involves powers of the sines of the angles, in a decreasing sequence from d-1 to zero. This means that the angles should be simulated according to their respective sine-power densities. However, except for the d=3 case, where simulating from the density sin(φ) is straightforward by inverse cdf, i.e. φ=acos(1-2u), the cdfs for the higher powers are combinations of sines and cosines, and as such are not easily inverted. Take for instance the eighth power:</p>
<p>F⁸(φ)=(840 φ – 672 sin(2 φ) + 168 sin(4 φ) – 32 sin(6 φ) + 3 sin(8 φ))/3072</p>
<p>While the densities are bounded by sin(φ), up to a constant, and hence an accept-reject can be easily derived, the efficiency decreases with the dimension according to the respective ratio of the <a rel="nofollow" href="https://en.wikipedia.org/wiki/Wallis%27_integrals" target="_blank">Wallis’ integrals</a>, unsurprisingly. A quick check for d=4 shows that the Normal simulation+projection-by-division-by-its-norm is faster.</p>
<p><a rel="nofollow" href="https://xianblog.files.wordpress.com/2018/03/xvaledarjct.jpg" target="_blank"><img data-orig-file="https://xianblog.files.wordpress.com/2018/03/xvaledarjct.jpg?w=450&amp;h=450&amp;fit=450%2C450&amp;resize=450%2C450" data-orig-size="480,480" data-medium-file="https://xianblog.files.wordpress.com/2018/03/xvaledarjct.jpg?w=450&amp;h=450&amp;fit=450%2C450&amp;resize=450%2C450?w=300" data-attachment-id="35591" data-image-meta="{" data-image-description="" src="https://xianblog.files.wordpress.com/2018/03/xvaledarjct.jpg?w=450&h=450&fit=450%2C450&resize=450%2C450" data-permalink="https://xianblog.wordpress.com/2018/03/08/uniform-on-the-sphere-or-not/xvaledarjct/" data-large-file="https://xianblog.files.wordpress.com/2018/03/xvaledarjct.jpg?w=450&amp;h=450&amp;fit=450%2C450&amp;resize=450%2C450?w=450" alt="" data-image-title="xvaledarjct" sizes="(max-width: 450px) 100vw, 450px" width="450" data-comments-opened="1" class="wp-image-35591" data-recalc-dims="1" height="450"></a>Puzzling a bit further about this while running, I wondered at the simultaneous simulations from sin(φ), sin(φ)², sin(φ)³, &amp;tc., but cannot see a faster way to recycle simulations from sin(φ). Points (φ,u) located in-between two adjacent power curves are acceptable simulations from the corresponding upper curve but they need be augmented by points (φ,u) under the lower curve to constitute a representative sample. In the end, this amounts to multiplying simulations from the highest power density as many times as there are powers. No gain in sight… Sigh!</p>
<p><a rel="nofollow" href="https://xianblog.files.wordpress.com/2018/03/densin3.jpg" target="_blank"><img data-orig-file="https://xianblog.files.wordpress.com/2018/03/densin3.jpg?w=450&amp;h=450&amp;fit=450%2C450&amp;resize=450%2C450" data-orig-size="480,480" data-medium-file="https://xianblog.files.wordpress.com/2018/03/densin3.jpg?w=450&amp;h=450&amp;fit=450%2C450&amp;resize=450%2C450?w=300" data-attachment-id="35597" data-image-meta="{" data-image-description="" src="https://xianblog.files.wordpress.com/2018/03/densin3.jpg?w=450&h=450&fit=450%2C450&resize=450%2C450" data-permalink="https://xianblog.wordpress.com/2018/03/08/uniform-on-the-sphere-or-not/densin3/" data-large-file="https://xianblog.files.wordpress.com/2018/03/densin3.jpg?w=450&amp;h=450&amp;fit=450%2C450&amp;resize=450%2C450?w=450" alt="" data-image-title="densin3" sizes="(max-width: 450px) 100vw, 450px" width="450" data-comments-opened="1" class="wp-image-35597" data-recalc-dims="1" height="450"></a>However, a few days later, while enjoying the sunset over Mont Blanc(!), I figured out that <a rel="nofollow" href="https://stats.stackexchange.com/a/331850/7224" target="_blank">there exists a direct and efficient way</a> to simulate from these powers of the sine function. Indeed, when looking at the density of cos(φ), it happens to be the signed root of a Beta(½,(d-1)/2), which avoids the accept-reject step. Presumably this is well-known, but I have not seen this proposal associated with the uniform distribution on the sphere.</p>
	<p>
<div>
<div>To <strong>leave a comment</strong> for the author, please follow the link and comment on their blog: <strong><a href="https://xianblog.wordpress.com/2018/03/08/uniform-on-the-sphere-or-not/"> R – Xi'an's Og</a></strong>.</div>
<hr>
<a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers.com</a> offers <strong><a rel="nofollow" href="https://feedburner.google.com/fb/a/mailverify?uri=RBloggers">daily e-mail updates</a></strong> about <a rel="nofollow" href="https://www.r-project.org/" title="The R Project for Statistical Computing">R</a> news and <a rel="nofollow" href="https://www.r-bloggers.com/search/tutorial" title="R tutorials">tutorials</a> on topics such as: <a rel="nofollow" href="https://www.r-bloggers.com/search/data%20science" title="Data science">Data science</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Big%20Data" title="Big Data">Big Data, <a rel="nofollow" href="https://www.r-users.com/" title="R jobs">R jobs</a>, visualization (<a rel="nofollow" href="https://www.r-bloggers.com/search/ggplot2" title="ggplot and ggplot2 tutorials">ggplot2</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/boxplot" title="Boxplots using lattice and ggplot2 tutorials">Boxplots</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/map" title="Maps and gis">maps</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/animation" title="Animation in R">animation</a>), programming (<a rel="nofollow" href="https://www.r-bloggers.com/search/RStudio" title="RStudio IDE for R">RStudio</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/sweave" title="Sweave and literate programming">Sweave</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/LaTeX" title="LaTeX in R">LaTeX</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/SQL" title="SQL and databases">SQL</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/eclipse" title="Eclipse IDE for R">Eclipse</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/git" title="git and github, Version Control System">git</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/hadoop" title="Large data in R using Hadoop">hadoop</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Web+Scraping" title="Web Scraping of google, facebook, yahoo, twitter and more using R">Web Scraping</a>) statistics (<a rel="nofollow" href="https://www.r-bloggers.com/search/regression" title="Regressions and ANOVA analysis tutorials">regression</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/PCA" title="principal component analysis tutorial">PCA</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/time+series" title="Time series">time series</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/trading" title="finance trading">trading</a>) and more...
<img alt="" width="1" src="http://feeds.feedburner.com/~r/RBloggers/~4/d1gYJj_azv4" height="1"></a></div>