text:
# I Probability involving n dice

Tags:
1. Apr 26, 2016

### goraemon

I'm studying probability and am currently stuck on this question:

Let's say we have n distinct dice, each of which is fair and 6-sided. If all of these dice are rolled, what is the probability that there is at least one pair that sums up to 7?

I interpreted the above as being equivalent to the following:

1 - (Probability that there is no pair that sums up to 7)

So if I were to consider just one pair of dice, then the probability that the pair adds up to 7 is 1/6, I think?

So Pr(one pair doesn't add up to 7) = 5/6.

But then I'm kind of stuck on how to proceed. Because there are lots of possible pairs amongst the n die, and some of these pairs overlap...for example, (die1, die2) is a pair, (die1, die3) is a pair, and so on. So I don't know how to account for these overlaps.

I tried breaking down the problem into a number of cases where there is no way for any pair to add up to 7:

(1) All of the dice show exactly one number.
(2) All of the dice show exactly two numbers which do not add up to 7 -- e.g. All the dice show either 3 or 6. Or all the dice show 2 or 4. And so on...
(3) All of the dice show exactly three numbers e.g. (1, 2, 4), no two of which can possibly add up to 7.
(4) If all of the dice show 4 or more numbers, then there MUST exist a pair that adds up to 7, so I don't consider any of these cases.

I suppose that I could add up the probabilities for all of the above cases, then I'd have the total probability that no two dice add up to 7? But then, how do I compute these?

In fact, is there a better / easier approach than the one I have thought up? Thanks.

2. Apr 26, 2016

### PeroK

I might do it the way you suggest. I'd look at using induction on n.

3. Apr 26, 2016

### Staff: Mentor

Do you know a function which can tell you exactly how many pairs there are?

4. May 1, 2016

### zinq

I definitely agree that the way to begin this problem is to find the probability that no pairs sum to 7, which is of course 1 - the probability you are ultimately seeking.

It might help organize the calculation if you specify what the sample space is: the set {1, 2, 3, 4, 5, 6} raised to the cartesian power of n, just a fancy way of saying the set of "all ordered n-tuples consisting of any of those 6 numbers". Then specify what the probability is for each point of this space, and what subset of the space you are finding the probability of. Finally, break that subset into disjoint pieces, each of which you can find the probability of.

You have virtually done this already, but it doesn't hurt to conceive of what you have done systematically, since that helps to see how various pieces that you have broken your target event into are disjoint from each other.

It will help if you figure out how many possible pairs of these numbers do not sum to 7. And then how many triples of these numbers involve no pair that sums to 7.

5. May 2, 2016

### micromass

Let's do an example:

(2) You know that two numbers are shown from the $n$ dice. We have the following choices, and all these choices have the same probabilities:
a) 1 and 2
b) 2 and 3
c) 1 and 3
....
(How many choices do we have for two numbers to show?)

Let's analyse case (a). The other cases are similar. So we know that all dice show either $1$ or $2$. Since the dice are all distinct, this means there are $2^n$ ways of $n$ dice to show either $1$ or $2$. However, $2$ of these ways are to show only $1$ and $2$. So there are $2^n - 2$ ways to satisfy $(a)$.

6. May 5, 2016

### bpet

Let $p_k(n)$ be the probability that case k has occurred after n dice have been thrown. It is not hard to show that these probabilities satisfy a recurrence relation
$$\begin{pmatrix}p_1(n+1) \\ p_2(n+1) \\ p_3(n+1) \end{pmatrix} = \begin{pmatrix} \frac{1}{6} & 0 & 0 \\ \frac{2}{3} & \frac{1}{3} & 0 \\ 0 & \frac{1}{3} & \frac{1}{2} \end{pmatrix} \begin{pmatrix}p_1(n) \\ p_2(n) \\ p_3(n) \end{pmatrix}$$
with $p_1(1)=1, p_2(1)=0=p_3(1)$. Since the matrix is lower triangular and the eigenvalues (along the diagonal) are distinct, it's not overly difficult to find the matrix powers to solve the recurrence relation. After a bit of algebra we find that the probability of obtaining at least one pair adding to 7 after n throws of the dice is
$$1-(p_1(n)+p_2(n)+p_3(n)) = 1-\frac{4}{2^{n-1}}+\frac{4}{3^{n-1}}-\frac{1}{6^{n-1}}$$


text:
At what altitude do the effects of atmospheric drag completely disappear?

The spate of K-questions got me thinking: At what point does atmospheric drag disappear1 with the pressure from sunlight, solar wind, or other forces becoming overwhelmingly dominant? I'm guessing it's somewhere in the outer reaches of Earth's magnetosphere so the answer might not be a fixed altitude. So, at what point does Earth's atmosphere cease having a dominant effect on an object's orbit2? Put another way, at what point is drag from Earth's atmosphere no longer able to cause an object's orbit to lower?

1 Properly-speaking, all gases/molecules/etc. that are above the Earth's surface and that are gravitationally bound to the Earth are within Earth's atmosphere. Atmospheric drag can thus occur anywhere within Earth's gravitational well. I suggest though that at some point the atmosphere should be swept away by e.g. solar wind; if it is no longer be present it won't be able to induce drag.

2 I reason that at some point other forces will either cancel out said drag's effect on an object or apply a force that has a net effect of boosting its orbit. I am unsure about this (especially the second assertion), hence the question. (This is clearly wrong since external forces on opposite sides of the orbit will cancel out while drag forces are additive.)

• Looks like a good improvement overall. It subtly changed the question though so I may tweak it a bit. – Alex Hajnal Dec 28 '18 at 2:21
• Answer to the title "At what altitude do the effects of atmospheric drag completely disappear?" None, because it never completely disappears. – uhoh Dec 28 '18 at 4:46
• @uhoh I disagree. At one point one surely would get swept away by e.g. solar wind (albeit very slowly, maybe). At the very least, at some point the atmosphere itself should be swept away. At least that's my reasoning (which could well be flawed). – Alex Hajnal Dec 28 '18 at 4:53
• Drag does not disappear because density does not disappear. Even if there were other forces, they could be resisted by drag. I think you want to ask about the net effect of all the forces, but losses due to drag would still be there no matter what additional forces there are. It never "disappears" per se. – uhoh Dec 28 '18 at 4:57
• @uhoh I see what you're getting at. I'll update the question to clarify that. – Alex Hajnal Dec 28 '18 at 5:00

Helpful factoids about orbits. To raise an orbit, you need a force that "pushes from behind", who's net effect accelerates you in the prograde direction. For circular orbits you can call that the $$+\theta$$ or tangential or $$\mathbf{\hat{v}}$$ direction. Counterintuitively this push from behind slows you down by an amount equal to the delta-v imparted, but it does raise your orbit.

Drag actually speeds you up by pushing against you, but lowers your orbit.

Aerodynamic drag always pushes in the exact opposite direction of your motion, because that's the definition. Aerodynamic effects on spacecraft in directions perpendicular to your velocity (radial or out-of-plane) are called lift.

The problem with comparing aerodynamic drag to effects from the solar wind or photon pressure is that those particles and photons are moving so much faster than the spacecraft's motion around the Earth, or the spacecraft plus Earth's motion around the Sun. These effects could be prograde at one point in the orbit and retrograde a half an orbit later, and ultimately substantially average out in some cases.

So it's a bit of an apples-to-oranges comparison.

However, in the plot below, the author has tried to show accelerations caused by various effects in an order-of-magnitude sort of way. For example, at an altitude of about 1300 km the acceleration due to solar radiation pressure is about the same size as that due to drag during maximum solar activity. Drag at high altitude is dramatically modulated by heating and expansion of the atmosphere by particles from the Sun, so actually both are caused by the Sun even though one also involves the atmosphere.

Below is borrowed from my answer to the question The sorting of perturbational effects by the power (used here as well).

I found the following plot in the book Satellite Orbits; Models, Methods, Applications by Oliver Montenbruck and Eberhard Gill, Springer, 2000. The figure and description can also be found in google books. It's a low quality snapshot but it's hard to capture a dozen different dependencies over 20 orders of magnitude without showing the whole thing.

Here is the bit of text that discusses the figure in more detail:

The effect of various perturbations a s a function of geocentric satellite distance is illustrated in Figure 3.1. For the calculation of the influence of atmospheric drag on circular low-Earth satellite orbits, exospheric temperatures between 500K and 2000K (cf. Sect. 3.5) have been assumed. The area-to-mass ratio used in the computation of non-gravitational forces is 0.01 m2/kg. For specially designed geodetic satellites like LAGEOS, the corresponding value may be smaller by one or two orders of magnitude. The perturbations due to various Geopotential coefficients Jn,m and the lunisolar attraction have been calculated from rule-of-thumb formulas by Milani et al. (1987). For the purposes of comparison it is mentioned that a constant radial acceleration of 10-11 km/s2 changes the semi-major axis of a geostationary satellite by approximately 1 m.

Aside from the aforementioned forces, various minor perturbations are considered in Fig. 3.1 which produce accelerations in the order of 10-15 to 10-12 km/s2. The most are due to the radiation pressure, resulting from the sunlight reflected by the Earth (albedo), as well as relativistic effects and the solid Earth tides.

• Would you support the assertion that there is persistent atmosphere (albeit at vanishingly low density) throughout Earth's gravity well? – Alex Hajnal Dec 28 '18 at 5:51
• @AlexHajnal I think it is tricky semantically; atmosphere is constantly escaping to infinity but as it interacts with UV light and the solar wind it becomes harder and harder to call it atmosphere and easier to call it plasma. I don't think there's an easy one-liner answer. But there's no magic wall where drag absolutely ceases. – uhoh Dec 28 '18 at 5:53
• I would consider plasma to be atmosphere in the context of this question. Excellent answer, BTW. – Alex Hajnal Dec 28 '18 at 5:57
• @AlexHajnal Thanks. I would encourage you to explore that further with a new question! – uhoh Dec 28 '18 at 6:00
• @AlexHajnal btw I switched $+\theta$ to MathJax since I'm already using it for $\mathbf{\hat{v}}$. Spelling out the word "delta-v" rather than using $\Delta v$ is really common, it's not for lack of access to a greek $\Delta$ character, it seems to be a real stand-alone term. en.wikipedia.org/wiki/Delta-v – uhoh Dec 28 '18 at 6:02

format: plain


text:
# Cathodic Reactions in Corrosion

|

The kinetics of the cathodic reactions often determines the corrosion rate of a metal or alloy in an aqueous electrolyte. The primary cathodic reactions in aqueous corrosion are listed below.

##### Proton reduction reaction (acids)

$$2{H^ + } + 2{e^ - } \leftrightarrow {H_2}(g) \tag{1} \label{eq:1}$$

##### Water reduction reaction (bases)

$$2{H_2}O + 2{e^ - } \leftrightarrow {H_2}(g) + 2O{H^ - } \tag{2} \label{eq:2}$$

##### Oxygen reduction reaction (acids)

$${O_2}(sol.) + 4{H^ + } + 4{e^ - } \leftrightarrow 2{H_2}O \tag{3} \label{eq:3}$$

##### Oxygen reduction reaction (bases)

$${O_2}(sol.) + 2{H_2}O + 4{e^ - } \leftrightarrow 4O{H^ - } \tag{4} \label{eq:4}$$

##### Equivalnet reactions

Equations $\eqref{eq:2}$ is equivalent to $\eqref{eq:1}$ and $\eqref{eq:4}$ to $\eqref{eq:3}$, assuming dissociation of water to $H^{+}$ and $OH^{-}$ and subtracting $OH^{-}$ from both sides.

#### References

D.A. Jones, "The Technology and evaluation of corrosion," in Principles and Prevention of Corrosion, 2nd ed., (Upper Saddle River, NJ: Prentice-Hall, 1996), p.5-9. ISBN: 0-13-359993-0. Buy at Amazon


text:
Remark 86.21.4 (Variant for adic-star). Let $P$ be a property of morphisms of $\textit{WAdm}^{adic*}$. We say $P$ is a local property if axioms (1), (2), (3) of Situation 86.21.2 hold for morphisms of $\textit{WAdm}^{adic*}$. In exactly the same way we obtain a variant of Lemma 86.21.3 for morphisms between locally adic* formal algebraic spaces over $S$.

Comment #1961 by Brian Conrad on

In the parenthetical at the start, replace "adic-star" with "adic-$\ast$".

In your comment you can use Markdown and LaTeX style mathematics (enclose it like $\pi$). A preview option is available if you wish to see how it works out (just click on the eye in the toolbar).

Unfortunately JavaScript is disabled in your browser, so the comment preview function will not work.

All contributions are licensed under the GNU Free Documentation License.

In order to prevent bots from posting comments, we would like you to prove that you are human. You can do this by filling in the name of the current tag in the following input field. As a reminder, this is tag 0ANH. Beware of the difference between the letter 'O' and the digit '0'.


text:
1 file | 385.27 KB

# On Borel summability and analytic functionals

(2013) 43(3). p.895-903
Author
Organization
Abstract
We show that a formal power series has positive radius of convergence if and only if it is uniformly Borel summable over a circle with center at the origin. Consequently, we obtain that an entire function $f$ is of exponential type if and only if the formal power series $\sum_{n=0}^{\infty}f^{(n)}(0)z^{n}$ is uniformly Borel summable over a circle centered a the origin. We apply these results to obtain a characterization of those Silva tempered ultradistributions which are analytic functionals. We also use Borel summability to represent analytic functionals as Borel sums of their moment Taylor series over the Borel polygon.
Keywords
Analytic functionals, Entire functions of exponential type, Borel summability, Borel polygon, Silva tempered ultradistributions

• Estrada-Vindas On Borel summability and analytic functionals.pdf
• full text
• |
• open access
• |
• PDF
• |
• 385.27 KB

## Citation

Chicago
Estrada, Ricardo, and Jasson Vindas Diaz. 2013. “On Borel Summability and Analytic Functionals.” Rocky Mountain Journal of Mathematics 43 (3): 895–903.
APA
Estrada, R., & Vindas Diaz, J. (2013). On Borel summability and analytic functionals. ROCKY MOUNTAIN JOURNAL OF MATHEMATICS, 43(3), 895–903.
Vancouver
1.
Estrada R, Vindas Diaz J. On Borel summability and analytic functionals. ROCKY MOUNTAIN JOURNAL OF MATHEMATICS. 2013;43(3):895–903.
MLA
Estrada, Ricardo, and Jasson Vindas Diaz. “On Borel Summability and Analytic Functionals.” ROCKY MOUNTAIN JOURNAL OF MATHEMATICS 43.3 (2013): 895–903. Print.
@article{4115876,
abstract     = {We show that a formal power series has positive radius of convergence if and only if it is uniformly Borel summable over a circle with center at the origin. Consequently, we obtain that an entire function \$f\$ is of exponential type if and only if the formal power series \${\textbackslash}sum\_\{n=0\}\^{ }\{{\textbackslash}infty\}f\^{ }\{(n)\}(0)z\^{ }\{n\}\$ is uniformly Borel summable over a circle centered a the origin. We apply these results to obtain a characterization of those Silva tempered ultradistributions which are analytic functionals. We also use Borel summability to represent analytic functionals as Borel sums of their moment Taylor series over the Borel polygon.},
author       = {Estrada, Ricardo and Vindas Diaz, Jasson},
issn         = {0035-7596},
journal      = {ROCKY MOUNTAIN JOURNAL OF MATHEMATICS},
keyword      = {Analytic functionals,Entire functions of exponential type,Borel summability,Borel polygon,Silva tempered ultradistributions},
language     = {eng},
number       = {3},
pages        = {895--903},
title        = {On Borel summability and analytic functionals},
url          = {http://dx.doi.org/10.1216/RMJ-2013-43-3-895},
volume       = {43},
year         = {2013},
}

Altmetric
View in Altmetric
Web of Science
Times cited:


text:
The Book of Lemmas Proposition 15 - Maple Help

 Main Concept Let $\mathrm{AB}$ be the diameter of a circle, $\mathrm{AC}$ a side of an inscribed regular pentagon, $D$ the middle point of the arc $\mathrm{AC}$. Join $\mathrm{CD}$ and produce it to meet $\mathrm{BA}$ produced in $E$; join $\mathrm{AC}$, $\mathrm{DB}$ meeting in $F$, and draw $\mathrm{FM}$ perpendicular to $\mathrm{AB}$. Then, $\mathrm{EM}$ = (  ). At the beginning of the proof, you might need to note this fact:   The inscribed angle theorem states that an angle q inscribed in a circle is half of the central angle 2q, that subtends the same arc on the circle.   A consequence of this is: For four consecutive points  on a circle, .
  Proof: Let $O$ be the center of the circle and join $\mathrm{DA}$, $\mathrm{DM}$, $\mathrm{DO}$, $\mathrm{CB}$. Now:                       ,    and,         , hence,                  Further, the triangles $\mathrm{FCB}$, $\mathrm{FMB}$, are equal in all respects. Therefore, in the triangles $\mathrm{DCB}$, $\mathrm{DMB}$, the sides $\mathrm{CB}$, $\mathrm{MB}$, being equal and $\mathrm{BD}$ common, while the angles $\mathrm{CBD}$, $\mathrm{MBD}$ are equal,   However,          so that,                   , and,                        .   Therefore,                 . Now, in the triangle $\mathrm{DMO}$,                              ,                              , (because ). Therefore,               ; hence,                        .   Again:                              Therefore, in the triangles $\mathrm{EDA}$, $\mathrm{ODM}$,                                            ,                                            , and the sides $\mathrm{AD}$, $\mathrm{MD}$ are equal.   Hence, the triangles are equal in all respects and,                                                  , Therefore,                                 .


text:
MAT244--2018F > Quiz-7

Q7 TUT 0601

(1/1)

Victor Ivrii:
(a) Determine all critical points of the given system of equations.

(b) Find the corresponding linear system near each critical point.

(c) Find the eigenvalues of each linear system. What conclusions can you then draw about the nonlinear system?

(d)  Draw a phase portrait of the nonlinear system to confirm your conclusions, or to extend them in those cases where the linear system does not provide definite information about the nonlinear system.
\left\{\begin{aligned} &\frac{dx}{dt} = 1 - y,\\ &\frac{dy}{dt} = x^2 - y^2. \end{aligned}\right.

Bonus: Computer generated picture

Zhiya Lou:
Solution

Jerry Qinghui Yu:
Solve $1-y = 0\Rightarrow y=1$ and $x^2 - y^2 = (x+y)(x-y) = 0$, substituting $y=1$, we have 2 critical points $(-1, 1), (1,1)$

Computing the Jacobian yields
\begin{align}
J =
\begin{bmatrix}
F_x & F_y\\
G_x & G_y
\end{bmatrix}
=
\begin{bmatrix}
0 & -1\\
2x & -2y
\end{bmatrix}
\end{align}
Plugging in both critical points we have 2 linear systems, first
\begin{align}
x' =
\begin{bmatrix}
0 & -1\\
-2 & -2
\end{bmatrix}
x
\end{align}
solving for eigenvalues yields $\lambda_1 = -1 + \sqrt{3}\, \lambda_2 = -1 - \sqrt{3}$, since $\sqrt{3} > 1$, we have $\lambda_1 > 0 > \lambda_2$, so we conclude locally $(-1, 1)$ is a saddle.
Plug in $(1,1)$
\begin{align}
x' =
\begin{bmatrix}
0 & -1\\
2 & -2
\end{bmatrix}
x
\end{align}
solving for eigenvalues yields $\lambda_1 = -1 + i\, \lambda_2 = -1 - i$, a complex conjugate with negative real parts, we conclude locally $(1, 1)$ is an asymptotically stable spiral.

Mengfan Zhu:
For this question,
there are two critical points: (-1,1) (1,1)
At (-1,1), it's an unstable saddle point.
At (1,1), it's an asymptotically stable spiral point.


text:
# Fast Mean Estimation with Sub-Gaussian Rates

We propose an estimator for the mean of a random vector in R^d that can be computed in time O(n^4+n^2d) for n i.i.d. samples and that has error bounds matching the sub-Gaussian case. The only assumptions we make about the data distribution are that it has finite mean and covariance; in particular, we make no assumptions about higher-order moments. Like the polynomial time estimator introduced by Hopkins, 2018, which is based on the sum-of-squares hierarchy, our estimator achieves optimal statistical efficiency in this challenging setting, but it has a significantly faster runtime and a simpler analysis.

• 13 publications
• 34 publications
• 58 publications
09/19/2018

### Mean Estimation with Sub-Gaussian Rates in Polynomial Time

We study polynomial time algorithms for estimating the mean of a heavy-t...
08/13/2019

### A Fast Spectral Algorithm for Mean Estimation with Sub-Gaussian Rates

We study the algorithmic problem of estimating the mean of heavy-tailed ...
11/17/2020

### Optimal Sub-Gaussian Mean Estimation in ℝ

We revisit the problem of estimating the mean of a real-valued distribut...
07/12/2020

### A spectral algorithm for robust regression with subgaussian rates

We study a new linear up to quadratic time algorithm for linear regressi...
05/10/2019

### Selectivity Estimation with Deep Likelihood Models

Selectivity estimation has long been grounded in statistical tools for d...
06/02/2020

### Robust and efficient mean estimation: approach based on the properties of self-normalized sums

Let X be a random variable with unknown mean and finite variance. We pre...
02/12/2018

### Dimension-free PAC-Bayesian bounds for the estimation of the mean of a random vector

In this paper, we present a new estimator of the mean of a random vector...

## 1 Introduction

Estimating the mean of a population given a finite sample is arguably the most fundamental statistical estimation problem. Despite the broad applicability and the fundamental nature of this problem, an estimator achieving the optimal statistical rate has only been discovered recently. However the optimal computational complexity of such an estimator is not well-understood.

In this paper, we are interested in obtaining high confidence estimates of the mean in the simple setting where only the existence of the covariance of the distribution is assumed. That is, we would like to find the smallest such that given samples from a distribution with mean our estimator satisfies:

 P{∥^X−μ∥≥rδ}≤δ.

To understand the inherent statistical limit of this problem, let us consider the simplified setting where the covariance is the identity. The most natural estimator for the mean of the population is the sample mean

. From the Central Limit Theorem, the distribution of

satisfies , and assuming this conclusion holds for any allows an  satisfying

 rδ=O⎛⎝√dn+√log1/δn⎞⎠.

[Cat12] shows that this

is the optimal statistical performance achievable under such mild assumptions. However, the above confidence interval only holds true asymptotically when the number of samples goes to infinity or when the distribution is sub-Gaussian. For finite sample results with a heavy-tailed distribution, applying Chebyshev’s inequality to the empirical mean gives only

 rδ=Ω(√dnδ).

The above bound is weaker than the one obtained by the Central Limit Theorem in two ways, the dependence on the failure probability

is polynomial in instead of logarithmic and the term depending on is multiplied by the dimensionality as opposed to being part of a smaller additive term. Unfortunately, [Cat12] also shows the above result is tight. That is, for any , there exists a distribution for which the bound guaranteed by Chebyshev’s inequality is optimal.

The poor performance of the empirical mean is due to its sensitivity to large outliers that occur naturally as part of the sample. The median-of-means framework was devised as a means of circumventing such difficulties. It was independently developed in the one dimensional case by

[NY83, JVV86, AMS99] and was later extended to the multivariate case by [HS16, LO11, Min15]. As part of this framework, the samples are first divided into batches and the mean of the samples is computed within each batch to obtain estimates . Each of these has mean

and variance

. The empirical mean is simply the mean of these estimates, which is sensitive to outliers. The median-of-means estimator instead is the geometric median of the estimates, which has greater tolerance to outliers. The success of the median-of-means estimator is due to the fact that it relies on only a fraction of estimates being close to the mean as opposed to all the estimates being close. [Min15] shows this gives an improved value of as follows:

 rδ=O⎛⎝√dlog1/δn⎞⎠.

The confidence interval guaranteed by the median-of-means estimator is better than the one for the empirical mean by improving the dependence on , but it is still poorer than we might expect from the Central Limit Theorem. Subsequent work attempting to bridge this gap achieves better rates than those guaranteed by the median-of-means but with stronger assumptions on the data generating distribution111A rate of is achieved under a fourth moment assumption on the distribution. ([JLO17]). The question of whether it was statistically feasible to obtain confidence intervals of the form guaranteed by the Central Limit Theorem was finally resolved by [LM19]. They devised an improved estimator, based on the median-of-means framework, called the median-of-means tournament, which achieves CLT-like confidence intervals. While the median-of-means estimator relies on the concentration of the number of close to the mean in Euclidean norm, the median-of-means tournament relies on the fact that along every direction , the number of close to the projection of the mean concentrates. The freedom to choose a different set of for each direction allow one to obtain a much smaller confidence interval than the one for the median-of-means estimator. In subsequent work, following the PAC-Bayesian approach of [Cat12], [CG17a] proposed a soft-truncation based estimator which obtains CLT-like confidence intervals provided one has access to estimates of the trace and spectral norm of the covariance matrix.

However, it is not known whether the estimators from [LM19, CG17a] are computationally feasible, as there are no known polynomial time algorithms to compute them. In contrast, the median-of-means and empirical mean can be computed in nearly-linear time ([CLM16]). To alleviate this computational intractability, [CG17b] proposed an efficient polynomial time estimator which achieves optimal statistical performance up to second order terms, assuming the existence of higher order moments. The question of computational tractability was subsequently resolved by [Hop18], who showed that an algorithm based on a sum-of-squares relaxation of the median-of-means tournament estimator achieves the statistically optimal CLT-like confidence intervals. However, the runtime of this algorithm is exorbitantly large222Assuming standard runtimes of the Interior Point method for semidefinite programming ([Ali95]) ().

In this paper, we propose a new algorithm with a reduced runtime——and a significantly simpler analysis. Our algorithm is a descent-based method that iteratively improves an estimate of the mean. The main challenge of such an approach is to estimate the descent direction. To this end, we crucially leverage the structure of the solutions to semidefinite programming relaxations of polynomial optimization problems designed to test whether a estimate is close to the mean. Our main contributions are twofold; we first show how exact solutions to the polynomial optimization problem furnish suitable descent directions and that such descent directions can also be efficiently extracted from relaxations of such problems and secondly, we show that these descent directions can be used in a descent style algorithm for mean estimation. Our paper is organized as follow: in Section 2, we present our main result, then in Section 3, as a warm-up, we devise a descent style algorithm for the case where we are given exact solutions to the polynomial optimization problems mentioned previously and prove that this algorithm achieves optimal statistical efficiency. This sets the stage for Section 4, where we present our main algorithm based on semidefinite relaxations of the previously defined polynomial optimization problems, leading to computationally efficient sub-Gaussian mean estimation.

## 2 Main result

Formally, our main result333The constants are explicit but we believe sub-optimal. is as follows:

###### Theorem 1.

Let be i.i.d. random vectors with mean and covariance . Then Algorithm 1 instantiated with Algorithms 4 and 5 and run with inputs , target confidence , stepsize and number of iterations returns a vector satisfying:

 ∥x∗−μ∥≤max⎛⎝ϵ,480000⎛⎝√TrΣn+√∥Σ∥log1/δn⎞⎠⎞⎠,

with probability at least .

We can make the following comments:

• Our estimator is both statistically optimal and computationally efficient. It achieves sub-Gaussian performance under minimal conditions on the distribution, and its runtime is . See Section 4.2 for details.

• The dependence of the number of iterations, , on can be avoided by initializing the algorithm with the median-of-means estimate. In this case, we can instead use and obtain the same guarantees, avoiding any dependence on the knowledge of .

• The estimator depends on the confidence level . [DLLO16] propose an estimator which works for a whole range of but for a restricted class of distributions.

• Our result does not explicitly depend on the dimension and our algorithm can be extended to a Hilbert space by working within the finite dimensional subspace containing the data points.

## 3 Warm-up

We present in this section a simple descent based algorithm. This algorithm is computationally inefficient but achieves the same guarantees of Theorem 1 with a much simpler analysis which nevertheless illustrates the main ideas behind the algorithm and proof of Theorem 1.

### 3.1 Intuition

We provide some intuition for our procedure, which iteratively improves an estimate of the mean. We first consider the simpler problem of testing whether a given point is close to the mean. We draw our inspiration from the main technical insight of [LM19], who show that along any direction, most of the bucket means, , are close to the mean, . Thus, to test whether a point, , is far from the mean, it is sufficient to check whether there exists a direction where most of the are far away from along that direction. This is formally expressed in the following polynomial optimization problem:

 maxk∑i=1bi b2i=bi ∥v∥2=1 bi⟨v,Zi−x⟩≥b2ir∀i∈[k] (MTE)

This polynomial problem over the set of variables and is parameterized by , the current estimate and the bucket means . Its polynomial constraints are encoding the number of beyond a distance from when projected along a direction . Intuitively, this program tries to find a direction so as to maximize the number of beyond a distance from along that direction. Here, we know from ([LM19]) that for an appropriate choice of , along all directions , a large fraction of the are close to the mean. Formally, for all directions , (see Corollary 1 ). Therefore this optimization problem has a large value when is far from the mean and can be used to certify this.

Strikingly, the direction returned by the solution of the above problem also contains information about the location of the mean when is chosen appropriately, which enables improvement of the quality of the current estimate. As illustrated in Figure 1, the direction returned by this optimization problem is strongly correlated with the vector joining the current point to the mean .

Algorithm 2 Distance Estimation 1:  Input: Data Points , Current point 2:   3:  Return:
Algorithm 3 Gradient Estimation 1:  Input: Data Points , Current point 2:   = Distance Estimation 3:   4:  Return:

Therefore, moving a small distance along the vector should intuitively take us closer to the mean. Given solutions to the polynomial optimization problem MTE, we may iteratively improve our estimate until no further change is necessary.

### 3.2 Algorithm

In this section we put the intuition provided previously into practice and propose a procedure that estimates the mean in the ideal situation where MTE can be exactly solved (the method is formally described in Algorithm 1):

1. First, following the median of means framework, the samples are divided into buckets and the mean of the samples within each bucket is computed as .

2. Second, the estimate of the mean is iteratively updated using a descent approach, based on the solution of MTE. As mentioned in Section 3.1, we need to run MTE with an appropriate choice of for the solution to be correlated with the direction . In the Distance Estimation step of our algorithm, we estimate a suitable choice of (see Algorithm 2). This value of is subsequently used in the Gradient Estimation step, to obtain an appropriate descent direction (see Algorithm 3).

From this point on, we refer to the solution of polynomial equations MTE as .

### 3.3 Analysis warm-up

In this simplified setting, we provide an analysis of our method and show that it obtains the same guarantees as those presented in Theorem 1. This is formally expressed in the following theorem for Algorithm 1 instantiated with Algorithms 2 and 3.

###### Theorem 2.

Let be i.i.d. random vectors with mean and covariance . Then Algorithm 1 instantiated with Algorithms 2 and 3 and run with inputs , target confidence , stepsize and number of iterations returns a vector satisfying:

 ∥x∗−μ∥≤max⎛⎝ϵ,108000⎛⎝√TrΣn+√∥Σ∥log1/δn⎞⎠⎞⎠,

with probability at least .

The main steps involved in the proof are the following:

1. Distance Estimation: We show that the Distance Estimation step in Algorithm 2 provides an accurate estimate of the distance of the current point from the mean. See Lemma 1.

2. Gradient Estimation: Next, we show that when is far away from the mean , the vector obtained by solving  MTE in Algorithm 3 is well aligned with the vector joining the current point to the mean . See Lemma 2.

3. Gradient Descent: Combining the previous two steps, we prove that we eventually converge to a good approximation to the mean.

In the proofs of our lemmas relating to the correctness of the Distance Estimation and the Gradient Estimation steps, we make use of the following assumption:

###### Assumption 1.

For the bucket means, , we have:

 ∀v∈Rd,∥v∥=1⇒∣∣{i:⟨Zi−μ,v⟩≥300(√TrΣ/n+√k∥Σ∥/n)}∣∣≤0.05k

The assumption is a formalization of the insight of ([LM19]), which shows that along all directions, , most of the bucket means are within a small radius of the true mean, , with high probability444This will be made precise in Corollary 1..

First, we prove that the Distance Estimation step defined in Algorithm 2 is correct.

###### Lemma 1.

Under Assumption 1, for all in the running of Algorihm 1, satisfies:

 ∣∣dt−∥xt−μ∥∣∣≤300(√TrΣ/n+√∥Σ∥k/n).
###### Proof.

Let . We first prove the lower bound . We may assume that , as the alternate case is trivially true. For , we can simply pick the vector where is the unit vector in the direction of . Under Assumption 1, we have that for at least points:

 ⟨Zi−xt,v⟩=⟨Zi−μ,v⟩+⟨μ−xt,v⟩≥∥xt−μ∥−r∗=r.

This implies the lower bound holds in the case where .

For the upper bound , suppose, for the sake of contradiction, there is a value of for which the optimal value of is greater than . Let be the solution of . This means that for of the , we have:

 ⟨Zi−μ,v⟩=⟨Zi−xt,v⟩+⟨xt−μ,v⟩≥r−∥xt−μ∥>r∗.

This contradicts Assumption 1 and proves the upper bound. ∎

Next, we prove the correctness of the Gradient Estimation step from Algorithm 3.

###### Lemma 2.

In the running of Algorithm 1, let us assume satisfies:

 ∥μ−xt∥≥1200(√TrΣ/n+√∥Σ∥k/n), (1)

and let denote the unit vector in the direction of . Then, under Assumption 1, we have that:

 ⟨gt,Δ⟩≥12.
###### Proof.

Let . We have, from the definition of , that for of the , . We also have, under Assumption 1, that for of the . From the pigeonhole principle, there exists a which satisfies both those inequalities. Therefore, for that , the lower bound from Lemma 1 implies

 ∥μ−xt∥−r∗≤dt≤⟨Zj−xt,gt⟩=⟨Zj−μ,gt⟩+⟨μ−xt,gt⟩≤r∗+∥μ−xt∥⟨Δ,gt⟩.

By rearranging the above inequality and using the assumption on in Eq. (1), we get the required conclusion. ∎

To control the probability that Assumption 1 holds, we assume the correctness of the following corollary of Lemma 7, formalizing the insight of ([LM19]):

###### Corollary 1.

Let be i.i.d. random vectors with mean and covariance . Furthermore, assume . Then we have for all such that :

 ∣∣{i:⟨Yi−μ,v⟩≥300(√TrΛ/k+√∥Λ∥)}∣∣≤0.05k

with probability at least .

By instantiating Corollary 1 with the , we see that Assumption 1 holds with high probability.

Finally, we put the results of Lemma 1, Lemma 2 and Corollary 1 together to prove Theorem 2.

###### Proof of Theorem 2.

Assume first that Assumption 1 holds. Let . To start with, let us define the set . We prove the theorem in two cases:

1. Case 1: None of the iterates lie in . In this case, note that by Lemma 1 and the definition of , we have:

 34∥xt−μ∥≤dt≤54∥xt−μ∥. (2)

Moreover, we have by the definition of the update rule of in Algorithm 1:

 ∥xt+1−μ∥2 =∥xt−μ∥2+12dt⟨xt−μ,gt⟩+d2t16≤∥xt−μ∥2−dt∥xt−μ∥4+d2t16 ≤∥xt−μ∥2−316∥xt−μ∥2+25256∥xt−μ∥2≤2325∥xt−μ∥2,

where we have used Lemma 2 for the first inequality and the inequalities in Eq. (2) for the second inequality. By iteratively applying the above inequality, we get the conclusion of the theorem in this case.

2. Case 2: At least one of the iterates lies in . Therefore, we have from Lemma 1:

 dt≤1500(√TrΣ/n+√∥Σ∥k/n).

We also have at the completion of the algorithm, from another application of Lemma 1:

 ∥x∗−μ∥−300(√TrΣ/n+√∥Σ∥k/n)≤d∗≤dt≤1500(√TrΣ/n+√∥Σ∥k/n).

By re-arranging the above inequality, we get the desired result.

By Corollary 1, Assumption 1 holds with probability at least and therefore, the conclusions from Case 1 and Case 2 hold with probability . ∎

Bearing in mind that the polynomial optimization problem MTE is non-convex, we consider a convex relaxation in the following section.

## 4 Efficient Algorithm for Mean Estimation

In this section, we define a semi-definite programming relaxation of the polynomial optimization problem MTE. We then design new Distance Estimation and Gradient Estimation algorithms that use the tractable solutions to the relaxation instead of the original polynomial optimization problem. We then use these solutions to update our mean estimate along the same lines as those from Section 3, albeit with some added technical difficulty. Finally, we provide the analysis of the method and prove Theorem 1.

### 4.1 The Semi-Definite Relaxation of Mte

Here, we propose a semidefinite programming relaxation of MTE, a variant of the Threshold-SDP from ([Hop18]). We first define a semidefinite matrix symbolically

indexed by , the variables and and denote by the vector :

 maxk∑i=1X1,bi X1,bi=Xbi,bi X1,1=1 d∑j=1Xvj,vj=1 ⟨vbi,Zi−x⟩≥Xbi,bir ∀i∈[k] X≽0 (MT)

Similar to the polynomial optimization MTE, this optimization problem is also parameterized by a vector , and a matrix . We refer to solutions of this program as with denoting the optimal value and denoting the optimal solution.

The main contribution of our paper is in showing that the solutions to the relaxed optimization problem MTE can be used to improve the mean estimate similar to those of MT.

### 4.2 Algorithm

To efficiently estimate the mean, we instantiate Algorithm 1 to use solutions of MT instead of MTE. The new Distance Estimation and Gradient Estimation procedures are stated in Algorithms 4 and 5.

As opposed to the polynomial optimization problem, solutions to the relaxation may not necessarily return a single vector but rather a semidefinite matrix which corresponds to the relaxation of

. This matrix may not uniquely determine a direction of improvement. We, therefore, parse the solution to isolate a provably good direction of improvement and use this to iteratively improve our estimate. It is noteworthy that the singular value decomposition does not provide a sign direction. Thankfully the correct orientation is easily ascertained using the data points.

To analyze the runtime of Algorithm 1 with Algorithms 4 and  5, we first note that the semidefinite relaxation has variables. However, by projecting all the data down to a subspace containing the bucket means, we may effectively reduce the number of variables to with an time pre-processing step. Therefore, we are now left with variables. The runtime of interior point methods for solving semidefinite programs with variables and constraints is ([Ali95]). Furthermore, a single call of the Distance Estimation procedure can be efficiently implemented using rounds of binary search on the parameter . Therefore, the total cost of a single call to Algorithm 4 is . Similarly, the total cost of a call to Algorithm 5 is . Since the cost of each iteration is dominated by a single call of Algorithm 4 and 5, the total cost per iteration is . Since, we only run iterations, the total cost of the Algorithm 1 instantiated with Algorithms 4 and  5 is .

### 4.3 Analysis

We now prove Theorem 1. We follow the same lines as the proof of Theorem 2, but with the added technical difficulties arising from the use of the semi-definite relaxation.

1. Distance Estimation: We show that the Distance Estimation step in Algorithm 4 provides an accurate estimate of the distance of the current point from the mean. See Section 4.3.1.

2. Gradient Estimation: Next, we show that when is far away from the mean , the vector output by Algorithm 5 is well aligned with the vector joining the current point to the mean . See Section 4.3.2.

3. Gradient Descent: Combining the previous two steps, we prove that we eventually converge to a good approximation to the mean. See Section 4.3.3.

The following assumption is required to prove the correctness of the Distance Estimation and Gradient Estimation steps:

###### Assumption 2.

For the bucket means, , let denote the set of feasible solutions for . Then, we have for all ,

 maxX∈Srk∑i=1Xbi,bi≤k20.

The above assumption is a strengthening of Assumption 1 for the case where we use MT instead of MTE. We use the following fact at several points in the subsequent analysis:

###### Remark 1.

Note that Assumption 2 implies Assumption 1.

#### 4.3.1 Distance Estimation Step

In this subsection, we analyze the Distance Estimation step from Algorithm 4. We show that an accurate estimate of the distance of the current point from the mean can be found. We begin by stating a lemma that shows that a feasible solution for can be converted to a feasible solution for with a reduction in optimal value.

###### Lemma 3.

Let us assume Assumption 2. Let be a positive semi-definite matrix, symbolically indexed by and the variables and . Moreover, suppose that satisfies:

 X1,1=1,Xbi,bi=X1,bi,d∑j=1Xvj,vj=1,k∑i=1Xbi,bi≥0.9k.

Then, there is a set of at least indices such that for all :

 ⟨Zi−μ,vbi⟩

and a set of at least indices such that for all , we have .

###### Proof.

Let . We prove the lemma by contradition. Firstly, note that is infeasible for as the optimal value for is less than (Assumption 2). Note that the only constraints of that are violated by are constraints of the form:

 ⟨Zi−μ,vbi⟩

Now, let denote the set of indices for which the above inequality is violated. We can convert to a feasible solution for by setting to the rows and columns corresponding to the indices in . Let be the matrix obtained by the above operation. We have from Assumption 2:

 0.05k≥k∑i=1X′bi,bi=k∑i=1Xbi,bi−∑i∈TXbi,bi≥0.9k−|T|,

where the last inequality follows from the fact that . By rearranging the above inequality, we get the first claim of the lemma.

For the second claim, let denote the set of indices satisfying . We have:

 0.9k≤k∑j=1Xbj,bj=∑j∈RXbj,bj+∑j∉RXbj,bj≤|R|+0.85k−0.85|R|⟹k3≤|R|.

This establishes the second claim of the lemma. ∎

The following lemma shows that if the distance between the mean and a point is small then the estimate returned by Algorithm 4 is also small.

###### Lemma 4.

Suppose a point satisfies . Then, under Assumption 2, Algorithm 4 returns a value satisfying

 d′≤7500(√TrΣ/n+√k∥Σ∥/n).
###### Proof.

Let and . Suppose that the optimal value of is greater than and let its optimal solution be . Let and denote the two sets whose existence is guaranteed by Lemma 3. From, the cardinalities of and , we see that their intersection is not empty. For , we have:

 0.85r′≤⟨Zj−x,vbj⟩=⟨Zj−μ,vbj⟩+⟨μ−x,vbj⟩

where the first inequality follows from the fact that and the fact that is feasible for and the last inequality follows from the inclusion of in and Cauchy-Schwarz.

By plugging in the bounds on and , we get:

 ∥x−μ∥>6075(√TrΣ/n+√k∥Σ∥/n).

This contradicts the assumption on and concludes the proof of the lemma. ∎

The next lemma shows that the distance between the mean and a point can be accurately estimated as long as is sufficiently far from .

###### Lemma 5.

Suppose a point satisfies . Then, under Assumption 2, Algorithm 4 returns a value satisfying:

 0.95~d≤d′≤1.25~d.
###### Proof.

Let us define the direction to be the unit vector in the direction of . From Assumption 1 (which is implied by Assumption 2), the number of satisfying is less than . Therefore, we have that for at least points:

 ⟨Zi−x,−Δ⟩=⟨x−μ+μ−Zi,Δ⟩=∥x−μ∥−300(√TrΣ/n+√k∥Σ∥/n)≥0.95~d.

Along with the monotonicity555See Lemma 8 in Appendix A. of in , this implies the lower bound.

For the upper bound, we show that the optimal value of is less than . For the sake of contradiction, suppose that this optimal value is greater than . Let be a feasible solution of that achieves . Let and be the two sets whose existence is guaranteed by Lemma 3 and be an element in their intersection. We have for :

 0.85(1.25~d) ≤Xbj,bj1.25~d≤⟨Zj−x,vbj⟩=⟨Zj−μ,vbj⟩+⟨μ−x,vbj⟩

where the first inequality follows from the inclusion of in and the last inequality follows from the inclusion of in and Cauchy-Schwarz. By re-arranging the above inequality, we get:

 Xbj,bj>(1.0625~d−~d)(300(√TrΣ/n+√k∥Σ∥/n))−1>1,

which is a contradiction. Therefore, we get from the monotonicity of (see Lemma 8), that and this concludes the proof of the lemma. ∎

In this section, we analyze the Gradient Estimation step of the algorithm. We show that an approximate gradient can be found as long as the current point is not too close to the mean . The following lemma shows that we obtain a non-trivial estimate of the gradient in Algorithm 5.

###### Lemma 6.

Suppose a point satisfies and let be the unit vector along . Then under Assumption 2, Algorithm 5 returns a vector satisfying:

 ⟨g,Δ⟩≥115.
###### Proof.

In the running of Algorithm 5, let denote the solution of . We begin by factorizing the solution into with the rows of denoted by , and . We also define the matrix in . From the constraints in MT, we have:

 Xbi,bi=∥ubi∥2≤1⟹∥ubi∥≤1,d∑j=1Xvj,vj=d∑j=1∥uvj∥2=∥Uv∥2F=1⟹∥Uv∥F=1.

Let and denote the sets defined in Lemma 3. Let . By noting that , we have for :

 0.85d∗≤⟨Zj−μ,vbj⟩+⟨μ−x,vbj⟩≤Xbj,bj300(√TrΣ/n+√k∥Σ∥/n)+u⊤bjUv(μ−x),

where the first inequality follows from the inclusion of in and the second from its inclusion in . We get by rearranging the above equation and using our bound on from Lemma 5:

 0.80∥μ−x∥≤0.85d∗≤Xbj,bj300(√TrΣ/n+√k∥Σ∥/n)+u⊤bjUv(μ−x). (3)

By rearranging Eq. (3), using Cauchy-Schwarz, and the assumption on :

 ∥Uv(μ−x)∥≥u⊤bjUv(μ−x)≥0.75∥μ−x∥.

We finally get that:

 ∥UvΔ∥≥0.75.

Now, we have:

Let be the top singular vector of . Note that and is also the top right singular vector of . We have that:

 0.75≤∥Uvy∥≤∥UvPΔy∥+∥UvP⊥Δy∥≤∥PΔy∥+∥UvP⊥Δ∥F≤∥PΔy∥+0.67.

This means that we have:

 |⟨y,Δ⟩|≥115.

Note that the algorithm returns either or . Firstly, consider the case where . From Assumption 1 (implied by Assumption 2), we have for at least points:

 ⟨Zi−μ,y⟩≤300(√TrΣ/n+√k∥Σ∥/n).

Therefore, we have for points:

 ⟨Zi−x,y⟩ =⟨Zi−μ,y⟩+⟨μ−x,y⟩ ≥−300(√TrΣ/n+√k∥Σ∥/n)+6000(√TrΣ/n+√k∥Σ∥/n)15>0.

This means that in the case where , we return which satisfies . This implies the lemma in this case. The case where is similar with used instead of . This concludes the proof of the lemma. ∎

The following lemma guarantees that Assumption 2 holds with high probability and is used analogously to Corollary 1 in the proof of Theorem 2:

###### Lemma 7.

Let be i.i.d. random vectors with mean and covariance and let denote the set of feasible solutions of . Then, we have for and :

 maxX∈Sk∑i=1Xbi,bi≤k20,

with probability at least .

The proof of the lemma is an application of standard empirical process theory and concentration inequalities ([LM19, Hop18]) and is proven in Appendix B.

The rest of the proof of Theorem 1 follows the same lines as that of Theorem 2 and is postponed to Appendix C.

## 5 Conclusion

In this paper, we proposed a computationally efficient estimator for the mean of a random vector which obtains the statistically optimal performance. This estimator has a significantly faster runtime together with a simpler analysis than previous works. Our algorithm is based on a descent method, where a current estimate of the mean is iteratively improved.

Considering the extension to M-estimation procedures ([BJL15, HS16, LM17]

) is a promising direction for further research, with as first step, the particular example of linear regression with heavy tailed noise and covariates (

[AC11]).

## References

• [AC11] J.-Y. Audibert and O. Catoni. Robust linear least squares regression. Ann. Statist., 39(5):2766–2794, 10 2011.

Interior point methods in semidefinite programming with applications to combinatorial optimization.

SIAM journal on Optimization, 5(1):13–51, 1995.
• [AMS99] N. Alon, Y. Matias, and M. Szegedy. The space complexity of approximating the frequency moments. Journal of Computer and System Sciences, 58(1):137–147, 1999.
• [BJL15] C. Brownlees, E. Joly, and G. Lugosi. Empirical risk minimization for heavy-tailed losses. Ann. Statist., 43(6):2507–2536, 12 2015.
• [BLM13] S. Boucheron, G. Lugosi, and P. Massart. Concentration Inequalities: A Nonasymptotic Theory of Independence. Oxford university press, 2013.
• [Cat12] O. Catoni. Challenging the empirical mean and empirical variance: a deviation study. Ann. Inst. Henri Poincaré Probab. Stat., 48(4):1148–1185, 2012.
• [CG17a] O. Catoni and I. Giulini. Dimension-free PAC-Bayesian bounds for matrices, vectors, and linear least squares regression. arXiv preprint arXiv:1712.02747, 2017.
• [CG17b] O. Catoni and I. Giulini. Dimension-free PAC-Bayesian bounds for the estimation of the mean of a random vector. NIPS 2017 workshop; (Almost) 50 shades of Bayesian learning: PAC-Bayesian trends and insights, 2017.
• [CLM16] M. B. Cohen, Y. T. Lee, G. Miller, J. Pachocki, and A. Sidford. Geometric median in nearly linear time. In

Proceedings of the Forty-eighth Annual ACM Symposium on Theory of Computing

, STOC ’16, 2016.
• [DLLO16] L. Devroye, M. Lerasle, G. Lugosi, and R. I. Oliveira. Sub-Gaussian mean estimators. Ann. Statist., 44(6):2695–2725, 2016.
• [Hop18] S. B Hopkins. Sub-Gaussian mean estimation in polynomial time. arXiv preprint arXiv:1809.07425, 2018.
• [HS16] D. Hsu and S. Sabato. Loss minimization and parameter estimation with heavy tails. J. Mach. Learn. Res., 17, 2016.
• [JLO17] E. Joly, G. Lugosi, and R. Oliveira. On the estimation of the mean of a random vector. Electron. J. Statist., 11(1):440–451, 2017.
• [JVV86] M. R. Jerrum, L. G. Valiant, and V. V. Vazirani.

Random generation of combinatorial structures from a uniform distribution.

Theoretical Computer Science, 43:169–188, 1986.
• [LM17] G. Lugosi and S. Mendelson. Risk minimization by median-of-means tournaments. Journal of the European Mathematical Society, to appear, 2017.
• [LM19] G. Lugosi and S. Mendelson. Sub-Gaussian estimators of the mean of a random vector. Ann. Statist., 47(2):783–794, 04 2019.
• [LO11] M. Lerasle and R. Oliveira. Robust empirical mean estimators. arXiv preprint arXiv:1112.3914, 2011.
• [LT91] M. Ledoux and M. Talagrand. Probability in Banach Spaces: Isoperimetry and Processes, volume 23. Springer Science & Business Media, 1991.
• [Min15] S. Minsker. Geometric median and robust estimation in Banach spaces. Bernoulli, 21(4):2308–2335, 2015.
• [Nes98] Y. Nesterov. Semidefinite relaxation and nonconvex quadratic optimization. Optimization methods and software, 9(1-3):141–160, 1998.
• [NY83] A. S. Nemirovsky and D. B. Yudin. Problem Complexity and Method Efficiency in Optimization. Wiley-Interscience Series in Discrete Mathematics. John Wiley & Sons, 1983.

## Appendix A Auxiliary lemma

###### Lemma 8.

For any and , the optimal value of is monotonically non-increasing in .

###### Proof.

The lemma follows trivially from the fact that a feasible solution of is also a feasible solution for for . ∎

## Appendix B Proof of Lemma 7

We first show that the optimal value of the semi-definite program MT satisfies a bounded-difference condition with respect to the ’s.

###### Lemma 9.

Let be any set of vectors in . Now, let be the same set of vectors with the vector replaced by . If and are the optimal values of and , we have:

 |m−m′|≤1
###### Proof.

Firstly, assume that is a feasible solution to . Now, let us define as:

 X′i,j={Xi,j if i,j≠bi0otherwise

That is is equal to except with the row and column corresponding to being set to . We see that forms a feasible solution to . Therefore, we have that:

 k∑j=1Xbj,bj=k∑j=1,j≠iX′bj,bj+Xbi,bi≤k∑j=1,j≠iX′bj,bj+1≤m′+1

where the bound follows from the fact that the sub-matrix of formed by the rows and columns indexed by and is positive semidefinite and the constraint that . Since the above series of equalities holds for all feasible solutions of , we get:

 m≤m′+1.

Through a similar argument, we also conclude that . Putting the above two inequalities together, we get the required conclusion. ∎

For the next few lemmas, we are concerned with the case where . Since we already know that the optimal SDP value satisfies the bounded differences condition, we need to verify that the expectation is small. As a first step towards this, we define the 2-to-1 norm of a matrix .

###### Definition 1.

The 2-to-1 norm of is defined as

 ∥M∥2→1=max∥v∥=1σi∈{±1}σ⊤Mv=max∥v∥=1∥Mv∥1

We consider the classical semidefinite programming relaxation of the 2-to-1 norm. To start with, we will define a matrix with the rows and columns indexed by and the elements and . The semidefinite programming relaxation is defined as follows:


text: So the fifth state, the Bose-Einstein Condensate is a bit quirky really. It is extremely hard to achieve such a state because you need to cool things down to a fraction of a degree above absolute zero (-273.15$^{0}C$). Getting things this cold is very difficult – even deep space has a temperature of 2.725$^{0}C$ which is not cold enough to achieve this temperature!


text:
# SL 360: Clear numeric entries on slide revisit

Hi there,

I have a few numeric entries on a slide and the slide property is set to 'Reset to initial state' on revisit. However, on revisit, the numeric entries still retain the entries from the previous attempt and are not resetting/clearing the fields. I tried resetting the numeric entry variable to value 0 on revisit, but this doesn't seem to fix it. Any help would be appreciated.

Thanks,

Prathiba

###### 2 Replies

You would also need triggers to reset the value of each variable.

Hi, thanks for the response. I added triggers to reset the value of each variable, however, that displayed 0 on revisit. Instead, we'd like to be able to clear all the fields upon revisit. I've added the below  JS script that I saw from this forum, but this now shows NaN upon revisit instead of cleared fields. Any suggestions will be great. Thanks so much!

var player = GetPlayer();
var empty1 = "";
player.SetVar("NumericEntry55", empty1);
player.SetVar("NumericEntry56", empty1);
player.SetVar("NumericEntry57", empty1);
player.SetVar("NumericEntry58", empty1);
player.SetVar("NumericEntry59", empty1);


text:
# How to Get a Job That You Don’t Hate

First of all, if you have not seen the movie Office Space, stop reading this blog post and watch it. I will wait.

Welcome back.

I am often asked to give career advice. This is strange since I don’t think I ever had a career. I had jobs, some terrible, some pretty good ones. I started a company in 2010. I am starting another one now. But a career, never. OK, so with that out of the way …

Earlier this month I was invited to discuss job search strategies with students in the MA in Statistics program at Columbia University. After the discussion, I posted the following blurb on my Facebook page.

Talking about careers in data analysis and stats with students in the MA in Statistics program at Columbia. My key messages: 1) if you want to work for banks, make sure you know what you are getting into; 2) think of an interview as a two way street: they interview you, you interview them; 3) if you hate your job, quit (if you can) and don’t worry about what it would look like on your resume; 4) don’t apply online, get a referral, go to meet ups, etc.; 5) learn some Bayesian stats — you will be a better human and know more than most of your peers.

I thought it would be useful to people if I elaborated on these a bit so here it goes.

## If you want to work for banks, make sure you know what you are getting into

A lot of students in the MA in Stats program want to work for banks. I am not sure why that is but it must have something do with the geography and expectations of high earnings. Whatever the motivation, it is a good idea to know what you are getting into. Not everyone hates working for banks, but in my experience, technical people who end up working there are not very happy. I think they find that the culture does not agree with them very much. My advice is always to ask to speak with your potential future peers and ask them, the future peers, about three things they love and three things they hate about their work. You would be surprised what you will learn. Having said that, I have met people on the “business” side of banks that absolutely love it. Like with anything else do your research and make your decisions based on conditional probabilities, not population averages.

## Think of an interview as a two-way street: they interview you, you interview them

This should be obvious, but most people don’t do it. The thing to recognize is that there is an inherent risk asymmetry between you and your prospective employer. You are just one candidate or employee out of many. They can make a mistake with you and they would probably be ok, but you are about to commit several years of your life to them (in expectation) and so you should be the one doing the interviewing! Of course, the realities of the sparse labor market is such that usually, you need them more than they need you, and so the roles are flipped. This fact, although daunting, should not deter you.

You want to find out what it would be like spending most of your waking hours at a job you do not yet have. This is not easy. To get started, make a simple two-category list: 1) culture; and 2) technical. For example, if you want flexible working hours, put that in the culture column and if you just must program in R, put that in the technical. Once you are done making the list, rank order the items. Do this before you take any interviews. After the interview, try to score the prospective employer along those dimensions. Where is the money column, you ask? That part is easy: know your minimum number and don’t be afraid to let them know what that is … but be reasonable, which means know what the market is paying and where you are on the skill / experience curve.

## If you hate your job, quit if you can and don’t worry about what it would look like on your resume

Some jobs are just plain awful. If you do what I recommended above, you will probably avoid most of those, but every now and again one will creep up on you. What to do? Quit! Sure, this is easier said than done, but at the very least immediately start looking for a new job and make some notes about how you were duped with this one. Introspection is a great tool and I use it often.

A friend of mine spent years working at a company for a horrible boss and even though he eventually quit he still has emotional anxiety over the whole affair. Life is way too short to work for assholes. Get out now. But what about the resume, you ask, and I answer: if you are a technical person, github (or something like that) is your resume.

## Don’t apply online, get a referral, go to meet-ups, and so on, but I am sorry I can’t refer you because I don’t know you

When I was working for a bank we had an opening for a business analyst. Now, here is the thing: business analyst does not analyze the business. What does she do? She writes requirements for a proposed piece of software. Anyway, that’s beside the point. When this job was posted by the HR department we received over 200 resumes! I don’t remember if we hired anyone, but you can imagine your chances of getting such a job. (Well, you can just compute them, but whatever.) The short story is, don’t apply online.

The best jobs I ever got were referred to me by my friends and classmates. Meetups are also some of the best places to get technical jobs. New York Statistical Programming Meetup is a great one for stats people and they often advertise jobs during their events. Another great way is to start contributing to some open source software. Where can you find great open source projects? Github, of course.

But Eric, why can’t you introduce me to some of those friends of yours that have all these great jobs? The truth is that they will not be my friends for much longer if I started doing that and you should not do it either. Your referral is a reflection on you — use it wisely and only introduce people you know well.

## Learn some Bayesian statistics — you will be a better human and know more than most of your peers

When I was getting my MA in Stats at CU, they did not have a masters level Bayesian class. This is a tragedy of modern statistical education, but things are getting better. My friend and co-founder Ben Goodrich is teaching an excellent Bayesian class for masters students in the QMSS program. The stats department also offers the class and Andrew Gelman teaches a PhD level Bayesian course. If you are not at Columbia, Coursera recently started offering Bayesian classes. This one looks particularly interesting.

So why all the hype about Bayes? It’s a long story, but here were my initial views on the subject. I now work exclusively in the Bayesian framework. In short, Bayes keeps me close to the reasons why I fell in love with statistics in the first place — it lets me model the world directly using probability distributions.

Even if you are not a statistician you should learn about conditional probabilities and Bayes rule so you do not commit many common fallacies such as the prosecutor’s fallacy, especially if you are a prosecutor.

## Bonus feature: why do you want a regular job anyway?

Recently I was on the Google hangouts call with a friend of mine who works as a contract programmer. His arrangement with  the company is that he works remotely. For most people remotely means working from home. For him, it means working from anywhere in the world. Right now he lives in a small apartment in Medellín, Columbia. He showed me the view from his window. It looks approximately like this:

To quote Tina Fey: I want to go to there.

The idea that an employer dictates both the hours during which you must work and the location of where the work must be performed is somewhat outdated. Sure, there are lots of jobs out there that legitimately require this kind of commitment, but it is no longer the norm. Take a look at that culture column I mentioned before and see where you stand relative to hours / location flexibility and choose accordingly.

## Note to people seeking H1B visa

A lot of people I speak with are in the US on F1 (student) visa. It is really tragic that the US does not award work visas to foreign graduates, but this is unlikely to change anytime soon. The common misconception is that you need to find a large company to sponsor your H1B (work visa). You do not. Lot’s of small companies can and do sponsor H1s. When I was working for a small startup in San Francisco in the mid-90s, we sponsored several H1Bs for Eastern European immigrants. The key is finding an experienced attorney who processes many applications and ask her for advice. Reputable attorneys will not charge you for the initial consultation.

# Talks, Lectures, and Workshops. What is the Difference?

I am about to go on a mini speaking tour and in preparation I am skimming Scott Berkun’s “Confession of a Public Speaker.” I like this book, but while reading it I realized that I will be giving two different types of “speeches”. Let’s call them talks and workshops, and even though in both cases the subject will be Stan, the audience’s expectations will be different and my presentation must reflect those differences. In particular, Scott’s book is a lot more relevant to talks than workshops.

Most inexperiences speakers assume that the people who come to their talks want to learn something and some people do have that expectation, but those are usually inexperienced consumers of talks. The truth is that it is very unlikely that you will learn something during a talk. Learning is a hard and active process and it is not going to happen by passively absorbing sound and light waves in a reclining position.  The most realistic goals for a talk is to inspire people to learn more about the subject. This is a difficult task for the presenter, but if you want to know how to do it well, I highly recommend Scott’s book.

A workshop is a different animal. As the name suggests, the participants will be working alongside the presenters and in so doing are hoping to come away with enough initial knowledge to jump start their own exploration. People who attend the workshop have already been inspired to learn more and the bar is therefore higher than during a talk. So what are the important attributes of a good workshop?

To think about that, image you are taking a technical class at a University. You are listening to a lecture. Are you at a talk or at a workshop? The listening part gives it away. Most likely you are at an uninspiring talk that should instead be a workshop. In order for the workshop to go well, here is my short of list of requirements:

1. Participants should have the required background at the right level of abstraction
2. If this is a computing workshop, participants already installed and tested the required software
3. Presenters have designed a series of exercises that gradually guide the audience through a set of hands on tasks each illuminating a different part of the subject
4. Participants have a chance to discuss the problem and their solutions with each other and with the instructors
5. There is a mechanism for the immediate feedback that tells the instructors if the majority have mastered the task

As a presenter, I can not control 1 and 2, but I must make it easy for people to assess their level of knowledge and software installation instructions must be clear.

Creating exercises is very time consuming, but I believe necessary for workshop style learning. Time for discussions can be weaved into the exercises and the output of the exercises can be shared with the rest of the class. Which brings me to the feedback mechanism, which is perhaps the most often overlooked aspect of the workshop.

I don’t have much experience with a feedback system during workshops, but I have used live surveys during talks and they work really well. For computing workshops, I would like to experiment with live code editors, where participants have a chance to post their code, their questions, and the error messages to the shared workspace. This would only work for moderate size groups, but I it seems to me that workshops should only be conducted in relatively small groups (say 50 people or fewer).

If you have any pointers on how to make the workshop experience better, feel free to post them in the comments.

# Thoreau: Thoughts on his Indictment and Defense

I never read Walden, not in its entirety anyway. I read most of the first chapter. It was dreadful. I still remember struggling to keep up with the narrative and wondering why is this such a big deal. Overall, I love the message of the simple life, civil disobedience, and living as one with nature. I do not love the apparently hypocritical obsession with seclusion and the disdain for all humanity. But this, of course, is a very shallow view of Thoreau. But then again, I do not have the patience to study him deeply. Fortunately, Kathryn Schultz and Jedediah Purdy do and offer an indictment of the man and somewhat halfhearted defense.

I really enjoyed reading both of these, but perhaps not surprisingly I found the indictment more convincing. The defence goes something like this. Sure, Thoreau was a hypocrite and an asshole, but we should not blame the message for the messenger (i.e. ad hominem or an opposite of blaming the messenger) even though in this case it happens to be the same person. I can get behind this argument. In science and in business there were and surely are lots of arrogant assholes, who nevertheless made important contributions. John Nash, despite a very favorable portrayal in the movie Beautiful Mind (the book is much less flattering), was not a very nice man. Steve Jobs was no sweetheart either. And so on. So, is Thoreau’s message important enough to stand on its own? That I am not qualified to answer, but a contrarian and anti-authoritarian in me wants to believe it that it is.

Thanks to Bryan Lewis for pointing me to these articles on his web page.

In 2012, I wrote a post about how to learn applied statistics without going to grad school. I still think that one does not have to spend a large amount of money to acquire the skills necessary for data analysis. What has changed for me personally is that am finding traditional statistical methods, call them classical or frequentist, or evolved classical based on the Stanford Statistical Learning school or whatever, somewhat unsatisfying.

These generally rely on maximum likelihood estimation (MLE) to generate point estimates and asymptotic properties of estimators to come up with confidence intervals. One of the main issues I have with this approach has nothing to do with MLE versus Bayesian full posterior per se. It has something to do with the fact that the Likelihood function is largely hidden from my view, although there are lots of others issues, some of which I hope to discuss when my understanding sufficiently progresses.  I am getting too comfortable just running glm(), ok, not glm() since there is no regularization there, but say glmnet or Random Forest or even bayesglm in R.  The latter is of course Bayesian, but still a black box.

I am not sure at this point if I am ready to abandon all the mathematical and algorithmic machinery of Lasso, Random Forests, Gradient Boosting Machines, and so on, but I would like to spend more time thinking and expressing models directly rather than running and tuning abstract algorithms. I am also quite certain I don’t want to write my own model fitting, sampling, and optimization procedures.

Since I would like to approach this problem in a Bayesian way, it also means that my goal is to get to the distribution of the parameter vector $\theta$ given data $y$, $p(\theta | y)$, the posterior. In the Bayesian framework, we still work with the likelihood function $p(y | \theta)$, but we are not trying to find some unique set of parameter values for which it is maximum (i.e. under which y are most likely.) Instead we want a complete picture of the uncertainty in our parameters that is supported by the data $y$, our choice of the model (i.e. likelihood, which as Andrew Gelman likes to point out is a form of prior knowledge) and knowledge about the parameters (prior distribution) without relying on asymptotic properties of estimators. In short:

Getting from prior to posterior is hard work unless they happen to be in the same family, which is rarely the case in the wild. The natural question then is where to start. Short of coding everything from scratch, which would be a very long project, even if I knew how to do it, two types of tools are in order: a probabilistic language capable of expressing models, parameters, priors, and their relationships and an MCMC sampler that can get us to the posterior distributions numerically. For a while, the best bet was some flavor of the BUGS language which uses Gibbs. But the state of the art has moved away from Gibbs sampling.  All the cool kids these days are playing with Stan which uses a more efficient, Hamiltonian MCMC with NUTS sampler and supports a broader set of models.

To get a jump start on Stan programming, I recently attended a class on Bayesian Inference with Stan taught by Andrew Gelman, Bob Carpenter, and Daniel Lee (thanks to Jared Lander for organizing the class.) I learned a lot and I hope to continue my exploration into Stan and Bayes.

* Thanks to Bob Carpenter for looking over the draft of this post and providing helpful comments.

# Getting There

After a 20 minute ride to the Jamaica Queens from Penn Station on the LIRR I am at the JFK.  LIRR is really the best way to travel to JFK from the city.  My flight leaves at 10:20 pm and the check-in procedure is quick an painless.  I noticed that Lufthansa clerk did not check if I had a Kazakhstan visa.

It is approximately 7 hours to Frankfurt, one hour layover and another 5.5 hours to Astana.  The flight was smooth and the airline food was surprisingly bearable.  On the flight to Astana three Russian woman were sitting next to me talking about a Government congress that they were either organizing or actively participating in; their enthusiasm for navigating the local bureaucracy was noticeable but alas not infectious. Their command of the English language was impressive.

We landed in Astana at 12:30 AM local time, right on schedule.  Anton, a grandson of repressed Polish people living in Russia during the Stalin era, met me at the gate.  20 minutes later we were in downtown Astana.

First day in Astana

A few hours of sleep and I was up and ready to soak in the local culture.  After acquiring a local sim card and telephone for $26 (4,000 ten) we visited the Pyramid and Court of Independence. Astana has an architectural plan, displayed on the third floor of the Castle of Independence spanning 18 years through 2030. The plan is to build a capital of Asia that would rival it seems the Arab Emirates. The Pyramid is a magnificent structure taking the shape according to its name. Several floors displaying (artwork?) culminate with a garden, leading to a small auditorium at the apex. There, the world religious congress is being held every three years and representatives from world religions meet in a friendly atmosphere. This was the first sign of religious tolerance and celebration of cultural diversity that I have encountered in Astana. It was one of many. After visiting the Pyramid Anton and I have taken a short walk across the street to visit the Court of Independence. It is another enormous structure of marble inside and out. The interior looks so polished that the floor ornaments reflect light making me feel like I am walking on an opaque mirror. An excursion organized by the staff takes us along the atrium to a set of stairs leading to the second floor. There, a small art gallery displays the works of local, American, and French painters, none familiar to me. I am particularly fond of the local painters depicting the nomadic lifestyle of the region. # Bayterek Bayterek was the first building that marked the beginning of the massive construction project undertaken by the president of Kazakhstan. Incidentally, his pictures and portraits are splattered all over Astana in museums, road signs, and government buildings. Normally this type of propaganda is a sign of deep political corruption, which may indeed be the case, but in Kazakhstan common people seem to have a lot of respect for the president. When I asked a tour guide at the top of Bayterek what she thought about the planned expansion of the city, she said that it would continue due to the will and ideas of the president, “may god bless his soul”. Her opinion was seconded by Anton, albeit not as enthusiastically. He thinks that even though Nazarbaev is the king of Kazakhstan, his policies of pluralism and his support of national diversity enables the peoples from 150 or so religions and nationalities to coexist in peace and mutual respect. This is not a small feat and the president should be applauded for this, irrespective of his motives. # Visit to the Chabad House The tour guide at Bayterek mentioned that among other religious institutions, Astana has a large synagogue. I was surprised to hear that, since I knew that very few Jews live in Kazakhstan. Naturally, I was excited to see it with my own eyes. The Astana Synagogue is in fact an impressive building of white and blue on the (Left) side of the city (Astana is divided in two informal parts by the river Ishum. The new construction is on the right bank of the river if you are facing (North?)) This Synagogue is fact part of Habad, a Hasidic organization whose mission is to attract more Jews into the religious practice. It is the only branch of orthodox Jewry that is welcoming of secular Jews. The building is guarded by a perimeter fence and a security guard at the front gate. After informing the guard that I would like to visit the Synagogue, he politely motions me inside and leads me up the stairs into the entrance. Anton follows me in. The Synagogue is empty at 2 PM in the afternoon. The guard informs the Rabbi of our arrival and in less than 30 seconds he appears smiling, sporting an obligatory long beard, white shirt, and a Hasidic hat. He is clearly happy to see us, and after inquiring if I were Jewish (he did not ask Anton, who looks very Slavic) hands both of us a yarmulke and invites us into the prayer hall. My guess is that the room sits approximately 50-100 people at the main level and another 20 or so on each balcony dividing the space into separate sections for men and women. Men can not focus on higher thoughts and prayers when women are too close, he explains. I think he is right about that. After he finished explaining the building attributes, I see that his excitement is starting to grow as he informs me that a big mitzvah is about to be performed. He runs to the back and promptly returns holding a black box containing a tefillin. I am a little bit squeamish, but he is so persistent that I submit and let him wrap my left hand with a flat leather rope, while placing a wooden box on my head. I repeat his words of prayer as he performs the procedure. When he is done, he informs me that he knows of 200-300 Jews that live in the city. For a city of 700,000 this is indeed a minority. Regular services usually gather 3-6 people and they only achieve minyan on Saturdays. Despite of a small congregation they publish a woman’s magazine and Jewish calendar. After giving Anton and I a copy of each, and asking us if we want something to eat or drink he politely escorted us back to the entrance. Please, come back whenever you are in town, he said to me. You as well, he said to Anton. # Trip to ALZHIR Before coming to Kazakhstan I knew that if I only could visit one site in Astana it would be ALZHIR. The camp lies 20 km East of the city in a village that was called Molinovka (meaning made from raspberries). It was part of Soviet Gulag system from 1937 until shorty after Stalin’s death in 1955. The system is commonly referred to as Gulag Archipelago due to Aleksandr Solzhenitsyn. ALZHIR’s significance is not that it housed political prisoners convicted under Article 58 of the NKVD code, but that it housed their wives and children under the age of 3 (children that were 3 years old or older were separated from their mothers and their siblings and sent to a ‘det dom’: a children’s house that is roughly equivalent to a penitentiary for juvenile delinquents). The road to ALZHIR (Akmolinki Lager’ Zhen Izmennikov Rodini: Akmolinsky Camp for the Wives of Traitors to the Motherland) is a deteriorating Soviet era asphalt road with enough bumps and cracks to ruin the suspension system of a modern tank. After about 30 minutes of bumping and grinding we arrive at the camp. As we approach the main building, the first item is in the front right – a Soviet rail car used in transporting women from all over the Soviet Union to ALZHIR. These wagons were manufactured in Odessa, Ukraine from 1928 to 1929 and women spent more than two months inside while being transported to ALZHIR. The wagon is approximately 7×3 meters (23×10 feet) and housed over 70 prisoners. Many had to stand and sleep on top of each other. To the left, a monument commemorating the victims. The barracks have been destroyed and the only surviving structure is the guard tower with the museum building in the back. Over 17,000 women were held in the camp, many have died from malnutrition and 18 hour work days. Their bodies were dumped into a mass grave behind the camp. Among the prisoners were noted Russian intellectuals and artists including Rachel Messerer, the mother of a famous Russian jewish ballerina Maya Mikhailovna Plisetskaya. Inside the museum the story of the prisoners unfolds. While still in Moscow, the women were told by the authorities that a rendezvous with their husbands was being arranged and asked to come to the prison building for the meeting. Many were so happy that they wore red dresses and makeup in anticipation. When they arrived, they were escorted to the ‘waiting area’. Many never saw freedom again. They were read their charges, asked to admit that their husbands were traitors, separated from their children, and packed into wagons for transport to ALZHIR. While at the camp they built their own barracks from straw and clay that could be found on the outskirts of the camp. The women also made their own clothes and later made uniforms for soldiers during the Second World War. One day, while collecting straw, a group of local Kazakh men and women came by. They watched them work and then started throwing what appeared to be small bricks and stones at the prisoners. The guards were laughing. Even the locals know that you are traitors, they said. One of the women picked up a small brick. It was made out of bread. When we left the camp, I asked Anton to stop in the middle of the steppe. In March the steppe is still covered with ice as temperatures can fluctuate between -30 (c) and -10 (c). I walked on the other side of the road and stepped onto the brittle snow. Even though it looked solid, my legs went in knee deep. For a while I just stood there breathing cold air and getting lost in the serenity of infinite spaces stretching for hundreds of miles in all directions. The Kazakh flatland is calming and magnetic. How was it, asked Anton when I got back in the car. It was magical, I said. Finally, I felt like I was in Kazakhstan. # First Two Weeks of Writing Jacki and I just submitted the first two chapter to our publisher, so I would like summarize early lessons learned (actually we submitted one chapter, but the editor decided to break the chapters in half; a decision that we fully support.) The chapters includes material on programming style (from R’s point of view), introduction to functions and functional programming, some information on S4 classes mostly from user’s perspective, vectorizing code, debugging and various methods of data access including web scraping and Twitter API. First the obvious. We underestimated the amount time required to produce the content. No surprises there. We spent too much time wrestling with the outline. Outlining seems to work well when I know my own writing style, but not so well otherwise. At some (earlier) point we should have just started writing and figured out the detailed chapter structure as we went along. I suspect this will change as we get deeper into the material, but only time will tell. What does need to be planned is the immediate section. For me it helps to have all the code written and all the visuals produced prior to starting writing. When I tried writing code on the fly, I struggled to make any meaningful progress. Lastly, it would have really helped if we read each other’s sections more carefully both in terms of synchronizing content and writing style. I hope that the final product does not read like the book was written by two people. Onto Chapter 2. # Getting Ready to Write a Book My co-author, Jacki Buros, and I have just signed a contract with Apress to write a book tentatively entitled “Predictive Analytics with R”, which will cover programming best practices, data munging, data exploration, and single and multi-level models with case studies in social media, healthcare, politics, marketing, and the stock market. Why does the world need another R book? We think there is a shortage of books that deal with the complete and programmer centric analysis of real, dirty, and sometimes unstructured data. Our target audience are people who have some familiarity with statistics, but do not have much experience with programming. Why did we not call the book Data Science blah, blah, blah…? Because Rachel and the Mathbabe already grabbed that title! (ok, kidding) The book is projected to be about 300 pages across 8 chapters. This is my first experience with writing a book and everything I heard about the process tells me that this is going to be a long and arduous endeavor lasting anywhere from 6 to 8 months. While undertaking a project of this size, I am sure there will be times when I will feel discouraged, overwhelmed, and emotionally and physically exhausted. What better vehicle for coping with these feelings than writing about them! (this is the last exclamation point in this post, promise.) So this is my first post of what I hope will become my personal diary detailing the writing process. Here is the summary of the events thus far. • A publisher contacted me on LinkedIn and asked if I wanted to write a book. • Jacki and I wrote a proposal describing our target market, competition, and sales estimates based on comparables. We developed an outline and detailed description of each section. • We submitted our proposal (to the original publisher and two other publishers) and received an approval to publish the book from Apress’ editorial board. (Apress was not the original publisher. More on that process after the book is complete.) We set up a tracking project on Trello (thanks Joel and the Trello team), created a task for every chapter, and a included a detailed checklist for each task. We have not completed all of the data analysis required for the book, so this is going to be an exercise in model building as well as in writing. If you have any advice about how to make the writing process better or if you think we are batshit crazy, please, post in the comments. I hope to write a book that we can be proud of. We have a great editorial team and a technical reviewer who is kind of a legend in the R/S world. They will remain anonymous for now, but their identities will be revealed as soon as they give me permission to do so. I am looking forward to learning about the writing process, about statistics, and about myself. Let the journey begin. # Joint Distributions or B Celebrity Sex in London I remember the first time the concept of joint probability distribution was introduced to me I found it completely unintuitive (like so many topics in probability), declared myself too stupid to get it, and considered giving up on statistics. The problem was that they used these silly gambling examples to demonstrate the concept. Flip of the coin this, roll of the die that. Urrg. Once I reset it in the context that I could relate too, everything became easier. So here we go. Suppose you are a B level celebrity in London. Then suppose that the probability of you having sex on any given day is ⅔ (it would be virtually 1, if you were an A level celebrity in LA, but that would not make for interesting example). Also, suppose that the probability of a cloudy day is ⅘. Here is how we write it in the language of probabilities. X is a random variable tracking our sex patterns. In this case, X can take on the values of X = Sex and X = NoSex. Y is a random variable tracking weather in London during the day. In our case Y = Cloudy or Y = Sunny. It should be obvious that P(X=Sex) or P(X=NoSex) = P(Y=Cloudy) or P(Y=Sunny) = 1. The probability of the whole sample space is 1 in both cases (you either have sex or not or it is either sunny or cloudy). So P(X=Sex) + P(X=NoSex) = ⅔ + ⅓ = 1. By the same token P(Y=Cloudy) + P(Y=Sunny) = ⅘ + ⅕ = 1. In this analysis, we assume that the sex is independent of the weather, an elusive assumptions and the one that should not be confused with correlation. I will try to make the differences clear when I discuss conditional independence. OK, here is the punch line. When we say joint probability distribution, we mean the probability of the combination of the events in question. For finite and discrete random variables such as the ones we are talking about, we can summarize the joint distribution using a table.  Sum = 1 P(X=Sex) = 2/3 P(X=NoSex) = 1/3 P(Y=Sunny) = 1/5 P(X=Sex,Y=Sunny) = 2/15 P(X=NoSex,Y=Sunny) = 1/15 P(Y=Cloudy) = 4/5 P(X=Sex,Y=Cloudy) = 8/15 P(X=NoSex,Y=Cloudy) = 4/15 In the above table when I write P(X=Sex, Y=Sunny) I mean the probability of both sex AND cloudy weather. (Do not confuse this with a notation P(X | Y), which means that we are only interested in Sex given a particular Weather outcome had already occurred.) This is why it is called joint probability. The probabilities listed in the margins of the table, are called … marginal. Each value in the cells is the product of two marginals. For instance the probability of having sex while it is sunny is ⅔ * ⅕ = 2/15. The cool thing is that if you observe only the joint probabilities you can easily calculate the marginals by summing across rows or columns (the reverse does not work – you can’t get joints from marginals, but you can get them from a good dealer in the Bronx.) So, if I observe you, the celebrity in question, for 150 days having sex during 20 sunny days and 80 cloudy days (you were not having sex during the other 50 days, sorry), I may conclude that the marginal probability of you having sex in London is 2/15 + 8/15 = 10/15 = ⅔ = P(X=Sex) = P(X=Sex,Y=Sunny)+ P(X=Sex,Y=Cloudy), which is in fact consistent with our assumptions. This can be formalized as follows: $P(X) = \sum_{Y}^{} P(X,Y)$ This rule is quite general (it is called the sum rule of probabilities) and it says that for any joint distribution X,Y to get back the probability of X we have to sum across all possible values of Y. Try summing across other rows columns and you will see that results are consistent there as well. # Updike’s Rabbit, Poincare, and the Art of Honest Writing Cover of Rabbit, Run As I am reading Rabbit, Run, I am slowly recognizing the literary genius of John Updike and I can not help but to draw parallels to the artists of the second kind — mathematicians. Updike does not use the tricks of literary construction that are so prevalent in the popular literature and modern blog writing. There is nothing wrong with clever literary construction of course. It makes the pages turn, it draws you in and leaves you asking for more. If you have read John Grisham’s Time to Kill (his first and best novel, I think), you know what I am talking about. The problem is that this kind of prose gets tiring after a while as you sort of feel like the author is consciously tricking you. Not so with Updike. His storyline is quite ordinary as are his characters. He does not leave you hanging at the end pages and paragraphs. He simply tells. The beauty of his writing, it seems to me, is that the prose itself is so cleverly nuanced, yet so vivid, that it infuses extraordinary qualities into ordinary events and actors. For example, from Rabbit, Run, describing a foreplay with a plump prostitute: As swiftly, he bends his face into a small forest smelling of spice, where he is out of all dimension, and where a tender entire woman seems an inch away, around a kind of corner. When he straightens up on his knees, kneeling as he is by the bed, Ruth under his eyes is an incredible continent, the pushed-up slip a north of snow. When reading Updike, the reading itself is an incredible experience, a total escape into the Updike dimension that is as insightful as it is unique. This kind of prose seems completely out of reach for mere mortals who need to resort to literary tricks. Cover via Amazon I get a similar feeling when reading Henri Poincare’s The Value of Science (in English translation) in that his understanding of mathematics is so deep that it feels almost untouchable, yet he simply tells without the drama of other popularizers of science like say Hawking (a brilliant man) or Mlodinow (also no slouch.) Not to be outdone by the literary types, Poincare’s narration is so beautiful that it makes me want to learn French just to read him in the original. Here is Poincare on the nuances of Number Theory: He is a savant indeed who will not take is as evident that every curve has a tangent; and in fact if we think of a curve and straight line as two narrow bands, we can always arrange them in such a way that they have a common part without intersecting And here he is again on the scientific motivation. The scientist does not study nature because it is useful; he studies it because he delights in it, and he delights in it because it is beautiful. If nature were not beautiful, it would not be worth knowing, and if nature would not be worth knowing, life would not be worth living. It was Poincare who noted that: A scientist worthy of his name, above all a mathematician, experiences in his work the same impression as an artist; his pleasure is as great and of the same nature. The curious intersection of art and science has been noted by many. The fact that science has its own aesthetic beauty is not a byproduct of the scientific method. As Poincare so eloquently points out, it is the reason for its existence. # To plot or to ggplot, that is not the question Producing informative and aesthetically pleasing quantitative visualizations is hard work. Any tool or library that helps me with this task is worth considering. Since I do most of my work in R, I have a choice of using plot, the default plotting library, a more powerful lattice package, and ggplot, which is based on the Grammar of Graphics. There is usually a tradeoff between the expressiveness of the grammer and the learning curve necessary to master it. I have recently invested 3 days of my life learning the ins and outs of ggplot and I have to say that it has been most rewarding. The fundamental difference between plot and ggplot is that in plot you manipulate graphical elements directly using predefined functions, whereas in ggplot you build the plot one layer at a time and can supply your own functions, although you can do quite a bit (but not everything) with a function called qplot, which abstracts the layering from the user and works similar to plot. And therefore qplot is exactly where you want to start when upgrading from plot. To demonstrate, the following R code partly visualizes the famous iris dataset containing Sepal and Petal measurements of three species of Iris flower using the built in plot function. par (mar=c(3,3,2,1), mgp=c(2,.7,0), tck=-.012, las=1) with(iris, plot(Sepal.Length, Sepal.Width, col=as.numeric(Species)+1, pch=20)) lbs = levels(iris$Species)
legend('topright', legend=lbs,
col=2:4, cex=0.7, pch=20, box.lwd=0.5, pt.cex=0.6)

One of the problems with plot is that the default plotting options are poorly chosen, so the first line of code fixed the margins, tick marks, and the orientation of the y axis tick labels.  The parameter col=as.numeric(Species) + 1 fixes the color offset at Red as opposed to the default Black.  Type palette() at the R prompt to see the default color vector.

The last complication is that plot does not draw the legend for you; it must be specified by hand.  And so, if you run the above code in R, you should get the following output.

It took a little bit of work, but the output looks pretty good.  Following is the equivalent task using ggplot’s qplot function.

qplot(Sepal.Length, Sepal.Width, data = iris, colour = Species, xlim=c(4,8))

As you can see, ggplot chooses a lot more sensible defaults and in this particular case, the interface for specifying the intent of the user is very simple and intuitive.

A final word of caution.  Just like a skier who sticks to blue and green slopes is in danger of never making it out of the intermediate hell, so is the qplot user will never truly master the grammar of graphics.  For those who dare to use a much more expressive ggplot(…) function, the rewards are well worth the effort.

Here are some of the ggplot references that I found valuable.