<p><div><div><div><div data-width="100" data-show-faces="false" data-send="true" data-layout="button_count" data-height="21" data-href="https://www.r-bloggers.com/lr03-residuals-and-rmse/"></div></div><div><a data-url="https://www.r-bloggers.com/lr03-residuals-and-rmse/" data-counturl="https://www.r-bloggers.com/lr03-residuals-and-rmse/" data-count="horizontal" data-text="LR03: Residuals and RMSE" class="twitter-share-button" data-via="rbloggers" href="https://twitter.com/share"></a></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="http://rbertolusso.github.io/posts/LR03-residuals-RMSE"> intubate &lt;||>, XBRL, bioPN, sbioPN, and stats with R</a></strong>, and kindly contributed to <a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers)</a> </div>
<p>This is post #3 on the subject of linear regression, using R for computational demonstrations and examples. We cover here <em>residuals</em> (or prediction errors) and the <em>RMSE</em> of the prediction line.<br>
The first post in the series is <a target="_blank" rel="nofollow" href="https://rbertolusso.github.io/posts/LR01-correlation">LR01: Correlation</a>.</p>
<hr>
<p></p>
<ul>
<li>Acknowledgments: organization is extracted from:
<ul>
<li><em>Freedman, Pisani, Purves, Statistics, 4th ed.</em>, probably<br>
the best book on statistical thinking (it maybe has a total of 4-5 formulas). It is referred here as FPP.</li>
<li>A lot of what is good is due to Professor Rudy Guerra.</li>
</ul>
</li>
</ul>
<hr>
<h3>Previous</h3>
<p><a target="_blank" rel="nofollow" href="https://rbertolusso.github.io/posts/LR01-correlation">LR01: Correlation</a></p>
<p><a target="_blank" rel="nofollow" href="https://rbertolusso.github.io/posts/LR02-SDline-GoA-regression">LR02: SD line, GoA, Regression</a></p>
<h3>What we know so far</h3>
<ul>
<li><em>Correlation coefficient</em>: measure of <strong>linear</strong> association,<br>
or clustering about the SD line.</li>
</ul>
<p></p>
<ul>
<li>If the scatterplot of both variables is a <em>football-shaped</em> cloud of points,<br>
those points cluster about the SD line, and the relationship between both variables<br>
can be summarized by:
<ul>
<li>average of x-values, SD of x-values ($s_x$),</li>
<li>average of y-values, SD of y-values ($s_y$).</li>
<li>
<p>the correlation coefficient r.</p>
</li>
<li>Reminder: we are (still) not making any probabilistic assumption, so we are not treating<br>
any variable as random, so everything stays in lower case. We are only saying (for now)<br>
that the correlation coefficient <em>makes sense</em>, as a summary “descriptive” statistic,<br>
for football-shaped<br>
cloud of points (even non-random ones), and it makes progressively less sense the more we<br>
depart from a football-shaped cloud of points. Following this reasoning, we could be using<br>
standard deviation version FPP (dividing by $n$), or $s$ (dividing by $n-1$) because<br>
we are <em>not</em> using them to estimate $\sigma$ (the population standard deviation) as that<br>
would imply making probabilistic assumptions (by the way, <em>both</em> version of the standard<br>
deviations, even with the right probabilistic assumptions, are <em>biased</em> estimators of<br>
the population standard deviation).</li>
</ul>
</li>
</ul>
<p><img title="plot of chunk SD-line" alt="plot of chunk SD-line" src="https://i1.wp.com/rbertolusso.github.io/figure/source/2016-11-25-LR03-residuals-RMSE/SD-line-1.png?w=456" data-recalc-dims="1"></p>
<ul>
<li>Graph of Averages (GoA): discrete function defined by $\text{Ave}(y|x)$, depending on<br>
the size of the $x$-intervals chosen (as happens with histograms):</li>
</ul>
<p><img title="plot of chunk GoA" alt="plot of chunk GoA" src="https://i0.wp.com/rbertolusso.github.io/figure/source/2016-11-25-LR03-residuals-RMSE/GoA-1.png?w=456" data-recalc-dims="1"></p>
<p>The GoA should <em>always</em> work, also for non football-shaped cloud of points.</p>
<ul>
<li>
<p>The regression (or regression <em>function</em>), linear or non-linear, of $y$<br>
(dependent variable) on $x$<br>
(independent variable) is a <em>continous</em> function that is a “smoother” of the GoA<br>
(that is a “discrete” function). It always work (but we may not know its<br>
mathematical form). The regression is the <em>continous</em> version of $\text{Ave}(y|x)$<br>
(GoA is the <em>discrete</em> version). <em>If</em> we are making probabilistic assumptions<br>
(by considering $Y$ a random variable with given characteristics that we will<br>
see later), we say that the regression function<br>
is the conditional expectaton of $Y$ given $X=x$ (or, in mathematical notation,<br>
$\text{E}(Y|X=x)$).</p>
</li>
<li>
<p>The regression <em>line</em> is a smoothed version of the GoA that is <em>correct</em> only<br>
<em>when</em> the GoA <em>is linear</em>.</p>
</li>
<li>
<p><em>When</em> the GoA is linear, a <em>best fit</em> regression <em>line</em>, that passes through the<br>
<em>point of averages</em> $(\bar{x}, \bar{y})$, is found as:</p>
</li>
</ul>
<p></p>
<p>where the slope $b$ is:</p>
<p></p>
<p>and the intercept $a$ is:</p>
<p></p>
<p><img title="plot of chunk GoA-reg-line" alt="plot of chunk GoA-reg-line" src="https://i1.wp.com/rbertolusso.github.io/figure/source/2016-11-25-LR03-residuals-RMSE/GoA-reg-line-1.png?w=456" data-recalc-dims="1"></p>
<ul>
<li>The regression line is an <em>attenuated</em> version of the SD line (that has slope<br>
= $\frac{s_y}{s_x}$)</li>
</ul>
<p><img title="plot of chunk unnamed-chunk-1" alt="plot of chunk unnamed-chunk-1" src="https://i0.wp.com/rbertolusso.github.io/figure/source/2016-11-25-LR03-residuals-RMSE/unnamed-chunk-1-1.png?w=456" data-recalc-dims="1"></p>
<ul>
<li>
<p>T or F: Regression is the line $y = a + bx$.</p>
</li>
<li>
<p>T or F: The SD line is the regression line.</p>
</li>
</ul>
<h3>Residuals, or Prediction errors</h3>
<ul>
<li>
<p>The regression method can be used to <em>predict</em> (guess) y from x<br>
(or x from y…). However, do you expect <em>actual</em> values to satisfy<br>
the predictions (guesses)?</p>
</li>
<li>
<p>They will differ, but, by how much? (What do you think?)</p>
</li>
<li>
<p>We will answer this in a little bit: we first need to consider<br>
the <em>prediction errors</em>, or <em>residuals</em>.</p>
</li>
</ul>
<p>The <em>distance</em> of a point above (+) or below (-) the regression line is:</p>
<p>error = actual – predicted.</p>
<ul>
<li>In the Pearson’s data (which was the name of the red line plotted?), showing<br>
only distances for a subset of points, it looks like:</li>
</ul>
<p><img title="plot of chunk residuals" alt="plot of chunk residuals" src="https://i0.wp.com/rbertolusso.github.io/figure/source/2016-11-25-LR03-residuals-RMSE/residuals-1.png?w=456" data-recalc-dims="1"></p>
<p>The prediction error is usually called <em>residual</em>. For consistency with<br>
Sheather’s book, we will name it $\hat{e}_i$ (you can also find it as<br>
$\hat{\epsilon}_i$, $r_i$, and perhaps other variants),<br>
the <em>predicted value</em> (or <em>fitted value</em>) is denoted by $\hat{y}_i$.<br>
The actual value is the y-coordinate of the point, $y_i$.</p>
<p>(Note: we are still not making any probabilistic assumptions, but <em>hats</em> are<br>
used to denote estimators, that are random variables. Even if usually random<br>
variables are denoted with upper case, residuals, or prediction errors, are usually<br>
denoted with lower case (the lower case also include errors $e_i$).<br>
When making probabilistic assumptions,<br>
$\hat{y}_i$ will become $\hat{Y}_i$.<br>
Sorry about the abuse of notation.)</p>
<p>Mathematically:</p>
<p></p>
<p>where</p>
<p></p>
<div>
<pre><code><span>## Predicted points
</span><span>yhat</span><span> </span><span>&lt;-</span><span> </span><span>meany</span><span> </span><span>-</span><span> </span><span>r</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>meanx</span><span> </span><span>+</span><span> </span><span>r</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>x</span><span>
</span><span>## Prediction errors (residuals): actual - predicted
</span><span>ehat</span><span> </span><span>&lt;-</span><span> </span><span>y</span><span> </span><span>-</span><span> </span><span>yhat</span><span>
</span><span>hist</span><span>(</span><span>ehat</span><span>,</span><span> </span><span>breaks</span><span> </span><span>=</span><span> </span><span>&quot;Scott&quot;</span><span>,</span><span> </span><span>probability</span><span> </span><span>=</span><span> </span><span>TRUE</span><span>)</span><span>
</span></code></pre>
</div>
<p><img title="plot of chunk hist-residuals" alt="plot of chunk hist-residuals" src="https://i1.wp.com/rbertolusso.github.io/figure/source/2016-11-25-LR03-residuals-RMSE/hist-residuals-1.png?w=456" data-recalc-dims="1"></p>
<ul>
<li>The smaller the distance of each point to the line, the better<br>
the <em>fit</em> to the regression line.</li>
</ul>
<h3>The RMSE of the prediction line</h3>
<p>RMSE stands for Root Mean Square Error. This implies:</p>
<ol>
<li>squaring the errors, or residuals:</li>
</ol>
<p></p>
<ol>
<li>Finding their mean:</li>
</ol>
<p></p>
<ol>
<li>Taking square root of the resulting mean:</li>
</ol>
<p></p>
<p>Using R:</p>
<div>
<pre><code><span>(</span><span>RMSE</span><span> </span><span>&lt;-</span><span> </span><span>sqrt</span><span>(</span><span>mean</span><span>(</span><span>ehat</span><span>^</span><span>2</span><span>)))</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 2.434294
</code></pre>
</div>
<ul>
<li>According to FPP, “the root mean square error for regression says how far typical<br>
points are above or below the regression line. The RMSE is to the regression line<br>
as the SD is to the average. For instance, if the scatter diagram is football-shaped,<br>
about 68% of the points on the scatter diagram will be within one RMSE of the regression<br>
line, about 95% of then will be within 2 RMSE of the regression line”.</li>
</ul>
<p><img title="plot of chunk RMSE-bands" alt="plot of chunk RMSE-bands" src="https://i2.wp.com/rbertolusso.github.io/figure/source/2016-11-25-LR03-residuals-RMSE/RMSE-bands-1.png?w=456" data-recalc-dims="1"></p>
<div>
<pre><code><span>## Proportion of values contained between 1 RMSE
</span><span>print</span><span>(</span><span>mean</span><span>(</span><span>abs</span><span>(</span><span>ehat</span><span>)</span><span> </span><span>&lt;</span><span> </span><span>RMSE</span><span>))</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 0.7133581
</code></pre>
</div>
<div>
<pre><code><span>## Proportion of values contained between 2 RMSEs
</span><span>print</span><span>(</span><span>mean</span><span>(</span><span>abs</span><span>(</span><span>ehat</span><span>)</span><span> </span><span>&lt;</span><span> </span><span>2</span><span> </span><span>*</span><span> </span><span>RMSE</span><span>))</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 0.9573284
</code></pre>
</div>
<ul>
<li>
<p>Pretty good, no?</p>
</li>
<li>
<p>FPP also present a relation between RMSE and SDy (the one calculated by dividing by $n$)</p>
</li>
</ul>
<p></p>
<p>Let’s get it with R:</p>
<div>
<pre><code><span>## FPP version of sample standard deviation
</span><span>SD</span><span> </span><span>&lt;-</span><span> </span><span>function</span><span>(</span><span>x</span><span>)</span><span> </span><span>sqrt</span><span>(</span><span>mean</span><span>((</span><span>x</span><span>-</span><span>mean</span><span>(</span><span>x</span><span>))</span><span>^</span><span>2</span><span>))</span><span>
</span><span>(</span><span>RMSE2</span><span> </span><span>&lt;-</span><span> </span><span>sqrt</span><span>(</span><span>1</span><span>-</span><span>r</span><span>^</span><span>2</span><span>)</span><span>*</span><span>SD</span><span>(</span><span>y</span><span>))</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 2.434294
</code></pre>
</div>
<div>
<pre><code><span>RMSE</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 2.434294
</code></pre>
</div>
<ul>
<li>The interpretation of RMSE is that it represents the <em>typical</em> size of the residuals.</li>
</ul>
<h3>Baby steps towards probabilistic assumptions</h3>
<ul>
<li>When the cloud of points is football-shaped the distribution of $y$-values<br>
<em>inside a strip</em> is <strong>Normal</strong> with mean approximately equal to $\text{Ave}(Y|X=x)$<br>
and standard deviation approximately equal to RMSE.</li>
</ul>
<p><img title="plot of chunk GoA-strips" alt="plot of chunk GoA-strips" src="https://i2.wp.com/rbertolusso.github.io/figure/source/2016-11-25-LR03-residuals-RMSE/GoA-strips-1.png?w=456" data-recalc-dims="1"></p>
<p><img title="plot of chunk hist-strips" alt="plot of chunk hist-strips" src="https://i1.wp.com/rbertolusso.github.io/figure/source/2016-11-25-LR03-residuals-RMSE/hist-strips-1.png?w=456" data-recalc-dims="1"></p>
<div>
<pre><code><span>mean</span><span>(</span><span>ystrip66</span><span>)</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 67.65625
</code></pre>
</div>
<div>
<pre><code><span>sd</span><span>(</span><span>ystrip66</span><span>)</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 2.350964
</code></pre>
</div>
<div>
<pre><code><span>mean</span><span>(</span><span>ystrip70</span><span>)</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 69.76845
</code></pre>
</div>
<div>
<pre><code><span>sd</span><span>(</span><span>ystrip70</span><span>)</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 2.48953
</code></pre>
</div>
<ul>
<li>
<p>If we <em>only</em> know about $y$, and not $x$, what do we use as an<br>
approximation to the <em>spread</em> (<em>variation</em>) of $y$?</p>
</li>
<li>
<p>If we also know about $x$, what do we use as an approximation<br>
to the spread of $y$ in a strip?</p>
</li>
<li>
<p>How is RMSE in relation to SDy? (smaller, equal, bigger)?</p>
</li>
</ul>
<p>In the case of sample data, as happened with SD, RMSE is not used nowadays because it is biased<br>
if used as an estimator of the standard deviation $\sigma$ of the residuals $e_i$ of the population<br>
(at some point we will calculate which is the bias of its square as an estimator<br>
of the variance $\sigma^2$ of $e_i$). Of course, this has importance only when we<br>
are dealing with samples <em>and</em> making probabilistic assumptions.</p>
<p>The version that is used, called <em>Residual standard error</em>,<br>
is <em>also</em> a biased estimator<br>
of $\sigma$, but its square (called <em>Mean Square Error of the residuals</em><br>
and indicated by $\text{MS}_\text{res}$), is an unbiased estimator of the variance $\sigma^2$ of $e_i$.</p>
<p>To start getting used to the functions <code>lm</code>, <code>summary</code> and <code>anova</code> that we will use heavily for simple and multiple regression (we will learn more about them some posts from now), let’s see<br>
where these values can be found.</p>
<div>
<pre><code><span>fit</span><span> </span><span>&lt;-</span><span> </span><span>lm</span><span>(</span><span>y</span><span> </span><span>~</span><span> </span><span>x</span><span>)</span><span>
</span><span>(</span><span>sfit</span><span> </span><span>&lt;-</span><span> </span><span>summary</span><span>(</span><span>fit</span><span>))</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## 
## Call:
## lm(formula = y ~ x)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.8772 -1.5144 -0.0079  1.6285  8.9685 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 33.88660    1.83235   18.49   &lt;2e-16 ***
## x            0.51409    0.02705   19.01   &lt;2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.437 on 1076 degrees of freedom
## Multiple R-squared:  0.2513,	Adjusted R-squared:  0.2506 
## F-statistic: 361.2 on 1 and 1076 DF,  p-value: &lt; 2.2e-16
</code></pre>
</div>
<p>The <em>Residual standard error</em> is easily spotted at the beginning of the third line starting<br>
from the bottom. It is rounded to three decimals. The actual value with more decimals (still rounded)<br>
is:</p>
<div>
<pre><code><span>sfit</span><span>$</span><span>sigma</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 2.436556
</code></pre>
</div>
<p>that is slightly different than:</p>
<div>
<pre><code><span>RMSE</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 2.434294
</code></pre>
</div>
<ul>
<li>How to “correct” RMSE to agree with Residual standard error?</li>
</ul>
<div>
<pre><code><span>RMSE</span><span> </span><span>*</span><span> </span><span>sqrt</span><span>(</span><span>n</span><span>/</span><span>(</span><span>n</span><span> </span><span>-</span><span> </span><span>2</span><span>))</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 2.436556
</code></pre>
</div>
<ul>
<li>
<p>Now we are in peace. But… why $n – 2$ ???</p>
</li>
<li>
<p>We said that the Residual standard error is a (biased) estimator of $\sigma$. Its<br>
square, $\text{MS}_\text{res}$, is an unbiased estimator of $\sigma^2$.</p>
</li>
</ul>
<div>
<pre><code><span>sfit</span><span>$</span><span>sigma</span><span>^</span><span>2</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 5.936804
</code></pre>
</div>
<p>Its value, using <code>anova</code>:</p>
<div>
<pre><code><span>(</span><span>afit</span><span> </span><span>&lt;-</span><span> </span><span>anova</span><span>(</span><span>fit</span><span>))</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## Analysis of Variance Table
## 
## Response: y
##             Df Sum Sq Mean Sq F value    Pr(>F)    
## x            1 2144.6 2144.58  361.23 &lt; 2.2e-16 ***
## Residuals 1076 6388.0    5.94                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
</code></pre>
</div>
<p>is found on the intersection of the row <em>Residuals</em> and the column <em>Mean Sq</em>. The<br>
less-rounded value is:</p>
<div>
<pre><code><span>afit</span><span>$</span><span>&quot;Mean Sq&quot;</span><span>[</span><span>2</span><span>]</span><span>
</span></code></pre>
</div>
<div>
<pre><code>## [1] 5.936804
</code></pre>
</div>
<h3>Appendix</h3>
<h4>Code to load data and initial calculations</h4>
<div>
<pre><code><span>library</span><span>(</span><span>UsingR</span><span>)</span><span>
</span><span>## Pearson's data
## Father and son data
</span><span>data</span><span>(</span><span>father.son</span><span>)</span><span>
</span><span>x</span><span> </span><span>&lt;-</span><span> </span><span>father.son</span><span>$</span><span>fheight</span><span>
</span><span>y</span><span> </span><span>&lt;-</span><span> </span><span>father.son</span><span>$</span><span>sheight</span><span>
</span><span>n</span><span> </span><span>&lt;-</span><span> </span><span>length</span><span>(</span><span>x</span><span>)</span><span>
</span><span>## Summary statistics
</span><span>meanx</span><span> </span><span>&lt;-</span><span> </span><span>mean</span><span>(</span><span>x</span><span>)</span><span>
</span><span>meany</span><span> </span><span>&lt;-</span><span> </span><span>mean</span><span>(</span><span>y</span><span>)</span><span>
</span><span>r</span><span> </span><span>&lt;-</span><span> </span><span>cor</span><span>(</span><span>x</span><span>,</span><span> </span><span>y</span><span>)</span><span>
</span></code></pre>
</div>
<h4>Code to produce the SD line plot</h4>
<div>
<pre><code><span>plot</span><span>(</span><span>x</span><span>,</span><span> </span><span>y</span><span>,</span><span>
</span><span>xlim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>),</span><span> </span><span>ylim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>),</span><span>
</span><span>xaxt</span><span> </span><span>=</span><span> </span><span>&quot;n&quot;</span><span>,</span><span> </span><span>yaxt</span><span> </span><span>=</span><span> </span><span>&quot;n&quot;</span><span>,</span><span> </span><span>xaxs</span><span> </span><span>=</span><span> </span><span>&quot;i&quot;</span><span>,</span><span> </span><span>yaxs</span><span> </span><span>=</span><span> </span><span>&quot;i&quot;</span><span>,</span><span>
</span><span>main</span><span> </span><span>=</span><span> </span><span>&quot;Pearson's data. SD line&quot;</span><span>,</span><span>
</span><span>xlab</span><span> </span><span>=</span><span> </span><span>&quot;Father height in inches&quot;</span><span>,</span><span> </span><span>ylab</span><span> </span><span>=</span><span> </span><span>&quot;Son height in inches&quot;</span><span>)</span><span>
</span><span>axp</span><span> </span><span>&lt;-</span><span> </span><span>seq</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>,</span><span> </span><span>by</span><span> </span><span>=</span><span> </span><span>2</span><span>)</span><span>
</span><span>axis</span><span>(</span><span>1</span><span>,</span><span> </span><span>at</span><span> </span><span>=</span><span> </span><span>axp</span><span>,</span><span> </span><span>labels</span><span> </span><span>=</span><span> </span><span>axp</span><span>)</span><span>
</span><span>axis</span><span>(</span><span>2</span><span>,</span><span> </span><span>at</span><span> </span><span>=</span><span> </span><span>axp</span><span>,</span><span> </span><span>labels</span><span> </span><span>=</span><span> </span><span>axp</span><span>)</span><span>
</span><span>## Spread of the cloud
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>meanx</span><span>+</span><span>c</span><span>(</span><span>-1</span><span>,</span><span>1</span><span>)</span><span>*</span><span>sd</span><span>(</span><span>x</span><span>),</span><span> </span><span>col</span><span>=</span><span>&quot;blue&quot;</span><span>,</span><span> </span><span>lty</span><span>=</span><span>3</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>h</span><span>=</span><span>meany</span><span>+</span><span>c</span><span>(</span><span>-1</span><span>,</span><span>1</span><span>)</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>),</span><span> </span><span>col</span><span>=</span><span>&quot;blue&quot;</span><span>,</span><span> </span><span>lty</span><span>=</span><span>3</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>meanx</span><span>+</span><span>c</span><span>(</span><span>-2</span><span>,</span><span>2</span><span>)</span><span>*</span><span>sd</span><span>(</span><span>x</span><span>),</span><span> </span><span>col</span><span>=</span><span>&quot;blue&quot;</span><span>,</span><span> </span><span>lty</span><span>=</span><span>3</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>h</span><span>=</span><span>meany</span><span>+</span><span>c</span><span>(</span><span>-2</span><span>,</span><span>2</span><span>)</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>),</span><span> </span><span>col</span><span>=</span><span>&quot;blue&quot;</span><span>,</span><span> </span><span>lty</span><span>=</span><span>3</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>meanx</span><span>+</span><span>c</span><span>(</span><span>-3</span><span>,</span><span>3</span><span>)</span><span>*</span><span>sd</span><span>(</span><span>x</span><span>),</span><span> </span><span>col</span><span>=</span><span>&quot;blue&quot;</span><span>,</span><span> </span><span>lty</span><span>=</span><span>3</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>h</span><span>=</span><span>meany</span><span>+</span><span>c</span><span>(</span><span>-3</span><span>,</span><span>3</span><span>)</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>),</span><span> </span><span>col</span><span>=</span><span>&quot;blue&quot;</span><span>,</span><span> </span><span>lty</span><span>=</span><span>3</span><span>)</span><span>
</span><span>## Point of averages (center of the cloud)
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>meanx</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;green&quot;</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>h</span><span>=</span><span>meany</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;green&quot;</span><span>)</span><span>
</span><span>## SD line using equation and sd
</span><span>abline</span><span>(</span><span>a</span><span> </span><span>=</span><span> </span><span>meany</span><span> </span><span>-</span><span> </span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>meanx</span><span>,</span><span>
</span><span>b</span><span> </span><span>=</span><span> </span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>),</span><span> </span><span>col</span><span> </span><span>=</span><span> </span><span>&quot;blue&quot;</span><span>,</span><span> </span><span>lwd</span><span> </span><span>=</span><span> </span><span>4</span><span>)</span><span>
</span></code></pre>
</div>
<h4>Code to produce the Graph of averages plot</h4>
<div>
<pre><code><span>plot</span><span>(</span><span>x</span><span>,</span><span> </span><span>y</span><span>,</span><span>
</span><span>xlim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>),</span><span> </span><span>ylim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>),</span><span> </span><span>col</span><span>=</span><span>&quot;lightgrey&quot;</span><span>,</span><span>
</span><span>xaxt</span><span> </span><span>=</span><span> </span><span>&quot;n&quot;</span><span>,</span><span> </span><span>yaxt</span><span> </span><span>=</span><span> </span><span>&quot;n&quot;</span><span>,</span><span> </span><span>xaxs</span><span> </span><span>=</span><span> </span><span>&quot;i&quot;</span><span>,</span><span> </span><span>yaxs</span><span> </span><span>=</span><span> </span><span>&quot;i&quot;</span><span>,</span><span>
</span><span>main</span><span> </span><span>=</span><span> </span><span>&quot;Pearson's data. GoA&quot;</span><span>,</span><span>
</span><span>xlab</span><span> </span><span>=</span><span> </span><span>&quot;Father height in inches&quot;</span><span>,</span><span> </span><span>ylab</span><span> </span><span>=</span><span> </span><span>&quot;Son height in inches&quot;</span><span>)</span><span>
</span><span>axp</span><span> </span><span>&lt;-</span><span> </span><span>seq</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>,</span><span> </span><span>by</span><span> </span><span>=</span><span> </span><span>2</span><span>)</span><span>
</span><span>axis</span><span>(</span><span>1</span><span>,</span><span> </span><span>at</span><span> </span><span>=</span><span> </span><span>axp</span><span>,</span><span> </span><span>labels</span><span> </span><span>=</span><span> </span><span>axp</span><span>)</span><span>
</span><span>axis</span><span>(</span><span>2</span><span>,</span><span> </span><span>at</span><span> </span><span>=</span><span> </span><span>axp</span><span>,</span><span> </span><span>labels</span><span> </span><span>=</span><span> </span><span>axp</span><span>)</span><span>
</span><span>## Point of averages (center of the cloud)
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>meanx</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;green&quot;</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>h</span><span>=</span><span>meany</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;green&quot;</span><span>)</span><span>
</span><span>fround</span><span> </span><span>&lt;-</span><span> </span><span>72</span><span>
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>fround</span><span>+</span><span>c</span><span>(</span><span>-.5</span><span>,</span><span>.5</span><span>),</span><span> </span><span>lty</span><span>=</span><span>3</span><span>)</span><span>
</span><span>with</span><span>(</span><span>subset</span><span>(</span><span>father.son</span><span>,</span><span> </span><span>round</span><span>(</span><span>fheight</span><span>,</span><span>0</span><span>)</span><span> </span><span>==</span><span> </span><span>fround</span><span>),</span><span>
</span><span>points</span><span>(</span><span>fheight</span><span>,</span><span> </span><span>sheight</span><span>,</span><span> </span><span>pch</span><span>=</span><span>16</span><span>,</span><span> </span><span>col</span><span>=</span><span>ifelse</span><span>(</span><span>sheight</span><span> </span><span>></span><span> </span><span>meany</span><span> </span><span>+</span><span> </span><span>(</span><span>fheight</span><span>-</span><span>meanx</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>*</span><span>r</span><span>,</span><span> </span><span>&quot;darkgreen&quot;</span><span>,</span><span> </span><span>&quot;blue&quot;</span><span>)))</span><span>
</span><span>fround</span><span> </span><span>&lt;-</span><span> </span><span>64</span><span>
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>fround</span><span>+</span><span>c</span><span>(</span><span>-.5</span><span>,</span><span>.5</span><span>),</span><span> </span><span>lty</span><span>=</span><span>3</span><span>)</span><span>
</span><span>with</span><span>(</span><span>subset</span><span>(</span><span>father.son</span><span>,</span><span> </span><span>round</span><span>(</span><span>fheight</span><span>,</span><span>0</span><span>)</span><span> </span><span>==</span><span> </span><span>fround</span><span>),</span><span>
</span><span>points</span><span>(</span><span>fheight</span><span>,</span><span> </span><span>sheight</span><span>,</span><span>
</span><span>pch</span><span>=</span><span>16</span><span>,</span><span> </span><span>col</span><span>=</span><span>ifelse</span><span>(</span><span>sheight</span><span> </span><span>></span><span> </span><span>meany</span><span> </span><span>+</span><span> </span><span>(</span><span>fheight</span><span>-</span><span>meanx</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>*</span><span>r</span><span>,</span><span> </span><span>&quot;darkgreen&quot;</span><span>,</span><span> </span><span>&quot;blue&quot;</span><span>)))</span><span>
</span><span>## Graph of averages.
</span><span>sgav</span><span> </span><span>&lt;-</span><span> </span><span>with</span><span>(</span><span>father.son</span><span>,</span><span> </span><span>tapply</span><span>(</span><span>sheight</span><span>,</span><span> </span><span>round</span><span>(</span><span>fheight</span><span>,</span><span>0</span><span>),</span><span> </span><span>mean</span><span>))</span><span>
</span><span>## sgavnum &lt;- with(father.son, tapply(sheight, round(fheight,0), length))
</span><span>points</span><span>(</span><span>as.numeric</span><span>(</span><span>names</span><span>(</span><span>sgav</span><span>)),</span><span> </span><span>sgav</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;red&quot;</span><span>,</span><span> </span><span>pch</span><span>=</span><span>16</span><span>)</span><span>
</span><span>## text(as.numeric(names(sgav)), sgav, sgavnum, pos=3)
</span></code></pre>
</div>
<h4>Code to produce the Graph of averages + regression line plot</h4>
<div>
<pre><code><span>plot</span><span>(</span><span>x</span><span>,</span><span> </span><span>y</span><span>,</span><span>
</span><span>xlim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>),</span><span> </span><span>ylim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>),</span><span> </span><span>col</span><span>=</span><span>&quot;lightgrey&quot;</span><span>,</span><span>
</span><span>xaxt</span><span> </span><span>=</span><span> </span><span>&quot;n&quot;</span><span>,</span><span> </span><span>yaxt</span><span> </span><span>=</span><span> </span><span>&quot;n&quot;</span><span>,</span><span> </span><span>xaxs</span><span> </span><span>=</span><span> </span><span>&quot;i&quot;</span><span>,</span><span> </span><span>yaxs</span><span> </span><span>=</span><span> </span><span>&quot;i&quot;</span><span>,</span><span>
</span><span>main</span><span> </span><span>=</span><span> </span><span>&quot;Pearson's data. GoA and regression line&quot;</span><span>,</span><span>
</span><span>xlab</span><span> </span><span>=</span><span> </span><span>&quot;Father height in inches&quot;</span><span>,</span><span> </span><span>ylab</span><span> </span><span>=</span><span> </span><span>&quot;Son height in inches&quot;</span><span>)</span><span>
</span><span>axp</span><span> </span><span>&lt;-</span><span> </span><span>seq</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>,</span><span> </span><span>by</span><span> </span><span>=</span><span> </span><span>2</span><span>)</span><span>
</span><span>axis</span><span>(</span><span>1</span><span>,</span><span> </span><span>at</span><span> </span><span>=</span><span> </span><span>axp</span><span>,</span><span> </span><span>labels</span><span> </span><span>=</span><span> </span><span>axp</span><span>)</span><span>
</span><span>axis</span><span>(</span><span>2</span><span>,</span><span> </span><span>at</span><span> </span><span>=</span><span> </span><span>axp</span><span>,</span><span> </span><span>labels</span><span> </span><span>=</span><span> </span><span>axp</span><span>)</span><span>
</span><span>## Point of averages (center of the cloud)
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>meanx</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;green&quot;</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>h</span><span>=</span><span>meany</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;green&quot;</span><span>)</span><span>
</span><span>fround</span><span> </span><span>&lt;-</span><span> </span><span>72</span><span>
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>fround</span><span>+</span><span>c</span><span>(</span><span>-.5</span><span>,</span><span>.5</span><span>),</span><span> </span><span>lty</span><span>=</span><span>3</span><span>)</span><span>
</span><span>with</span><span>(</span><span>subset</span><span>(</span><span>father.son</span><span>,</span><span> </span><span>round</span><span>(</span><span>fheight</span><span>,</span><span>0</span><span>)</span><span> </span><span>==</span><span> </span><span>fround</span><span>),</span><span>
</span><span>points</span><span>(</span><span>fheight</span><span>,</span><span> </span><span>sheight</span><span>,</span><span> </span><span>pch</span><span>=</span><span>16</span><span>,</span><span> </span><span>col</span><span>=</span><span>ifelse</span><span>(</span><span>sheight</span><span> </span><span>></span><span> </span><span>meany</span><span> </span><span>+</span><span> </span><span>(</span><span>fheight</span><span>-</span><span>meanx</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>*</span><span>r</span><span>,</span><span> </span><span>&quot;darkgreen&quot;</span><span>,</span><span> </span><span>&quot;blue&quot;</span><span>)))</span><span>
</span><span>fround</span><span> </span><span>&lt;-</span><span> </span><span>64</span><span>
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>fround</span><span>+</span><span>c</span><span>(</span><span>-.5</span><span>,</span><span>.5</span><span>),</span><span> </span><span>lty</span><span>=</span><span>3</span><span>)</span><span>
</span><span>with</span><span>(</span><span>subset</span><span>(</span><span>father.son</span><span>,</span><span> </span><span>round</span><span>(</span><span>fheight</span><span>,</span><span>0</span><span>)</span><span> </span><span>==</span><span> </span><span>fround</span><span>),</span><span>
</span><span>points</span><span>(</span><span>fheight</span><span>,</span><span> </span><span>sheight</span><span>,</span><span>
</span><span>pch</span><span>=</span><span>16</span><span>,</span><span> </span><span>col</span><span>=</span><span>ifelse</span><span>(</span><span>sheight</span><span> </span><span>></span><span> </span><span>meany</span><span> </span><span>+</span><span> </span><span>(</span><span>fheight</span><span>-</span><span>meanx</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>*</span><span>r</span><span>,</span><span> </span><span>&quot;darkgreen&quot;</span><span>,</span><span> </span><span>&quot;blue&quot;</span><span>)))</span><span>
</span><span>## Graph of averages.
</span><span>sgav</span><span> </span><span>&lt;-</span><span> </span><span>with</span><span>(</span><span>father.son</span><span>,</span><span> </span><span>tapply</span><span>(</span><span>sheight</span><span>,</span><span> </span><span>round</span><span>(</span><span>fheight</span><span>,</span><span>0</span><span>),</span><span> </span><span>mean</span><span>))</span><span>
</span><span>## sgavnum &lt;- with(father.son, tapply(sheight, round(fheight,0), length))
</span><span>points</span><span>(</span><span>as.numeric</span><span>(</span><span>names</span><span>(</span><span>sgav</span><span>)),</span><span> </span><span>sgav</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;red&quot;</span><span>,</span><span> </span><span>pch</span><span>=</span><span>16</span><span>)</span><span>
</span><span>## text(as.numeric(names(sgav)), sgav, sgavnum, pos=3)
</span><span>
</span><span>## Regression line
</span><span>abline</span><span>(</span><span>a</span><span>=</span><span>meany</span><span>-</span><span>r</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>meanx</span><span>,</span><span> </span><span>b</span><span>=</span><span>r</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>),</span><span> </span><span>lwd</span><span>=</span><span>4</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;red&quot;</span><span>)</span><span>
</span></code></pre>
</div>
<h4>Code to produce the Residual plot</h4>
<div>
<pre><code><span>plot</span><span>(</span><span>x</span><span>,</span><span> </span><span>y</span><span>,</span><span>
</span><span>xlim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>),</span><span> </span><span>ylim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>),</span><span> </span><span>col</span><span> </span><span>=</span><span> </span><span>&quot;lightgrey&quot;</span><span>,</span><span>
</span><span>xaxt</span><span> </span><span>=</span><span> </span><span>&quot;n&quot;</span><span>,</span><span> </span><span>yaxt</span><span> </span><span>=</span><span> </span><span>&quot;n&quot;</span><span>,</span><span> </span><span>xaxs</span><span> </span><span>=</span><span> </span><span>&quot;i&quot;</span><span>,</span><span> </span><span>yaxs</span><span> </span><span>=</span><span> </span><span>&quot;i&quot;</span><span>,</span><span>
</span><span>main</span><span> </span><span>=</span><span> </span><span>&quot;Prediction errors (or Residuals)&quot;</span><span>,</span><span>
</span><span>xlab</span><span> </span><span>=</span><span> </span><span>&quot;Father height in inches&quot;</span><span>,</span><span> </span><span>ylab</span><span> </span><span>=</span><span> </span><span>&quot;Son height in inches&quot;</span><span>)</span><span>
</span><span>axp</span><span> </span><span>&lt;-</span><span> </span><span>seq</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>,</span><span> </span><span>by</span><span> </span><span>=</span><span> </span><span>2</span><span>)</span><span>
</span><span>axis</span><span>(</span><span>1</span><span>,</span><span> </span><span>at</span><span> </span><span>=</span><span> </span><span>axp</span><span>,</span><span> </span><span>labels</span><span> </span><span>=</span><span> </span><span>axp</span><span>)</span><span>
</span><span>axis</span><span>(</span><span>2</span><span>,</span><span> </span><span>at</span><span> </span><span>=</span><span> </span><span>axp</span><span>,</span><span> </span><span>labels</span><span> </span><span>=</span><span> </span><span>axp</span><span>)</span><span>
</span><span>## Regression line
</span><span>abline</span><span>(</span><span>a</span><span> </span><span>=</span><span> </span><span>meany</span><span> </span><span>-</span><span> </span><span>r</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>meanx</span><span>,</span><span> </span><span>b</span><span> </span><span>=</span><span> </span><span>r</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>),</span><span> </span><span>lwd</span><span> </span><span>=</span><span> </span><span>4</span><span>,</span><span> </span><span>col</span><span> </span><span>=</span><span> </span><span>&quot;red&quot;</span><span>)</span><span>
</span><span>## for (fn in sample(nrow(father.son), 20)) {
</span><span>for</span><span> </span><span>(</span><span>i</span><span> </span><span>in</span><span> </span><span>c</span><span>(</span><span>39</span><span>,</span><span> </span><span>158</span><span>,</span><span> </span><span>204</span><span>,</span><span> </span><span>479</span><span>,</span><span> </span><span>686</span><span>,</span><span> </span><span>808</span><span>,</span><span> </span><span>844</span><span>,</span><span> </span><span>851</span><span>,</span><span> </span><span>852</span><span>,</span><span> </span><span>1070</span><span>))</span><span> </span><span>{</span><span>
</span><span>## Actual point
</span><span> </span><span>points</span><span>(</span><span>x</span><span>[</span><span>i</span><span>],</span><span> </span><span>y</span><span>[</span><span>i</span><span>])</span><span>
</span><span>## Predicted point yhat_i
</span><span> </span><span>yhat_i</span><span> </span><span>&lt;-</span><span> </span><span>r</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>x</span><span>[</span><span>i</span><span>]</span><span> </span><span>+</span><span> </span><span>meany</span><span>-</span><span>r</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>meanx</span><span>
</span><span>points</span><span>(</span><span>x</span><span>[</span><span>i</span><span>],</span><span> </span><span>yhat_i</span><span>,</span><span> </span><span>pch</span><span>=</span><span>16</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;red&quot;</span><span>)</span><span>
</span><span>lines</span><span>(</span><span>rep</span><span>(</span><span>x</span><span>[</span><span>i</span><span>],</span><span> </span><span>2</span><span>),</span><span> </span><span>c</span><span>(</span><span>y</span><span>[</span><span>i</span><span>],</span><span> </span><span>yhat_i</span><span>))</span><span>
</span><span>## Prediction error (or residual): actual - predicted
</span><span> </span><span>ehat_i</span><span> </span><span>&lt;-</span><span> </span><span>y</span><span>[</span><span>i</span><span>]</span><span> </span><span>-</span><span> </span><span>yhat_i</span><span>
</span><span>text</span><span>(</span><span>x</span><span>[</span><span>i</span><span>],</span><span> </span><span>(</span><span>y</span><span>[</span><span>i</span><span>]</span><span> </span><span>+</span><span> </span><span>yhat_i</span><span>)</span><span>/</span><span>2</span><span>,</span><span> </span><span>round</span><span>(</span><span>ehat_i</span><span>,</span><span> </span><span>2</span><span>),</span><span> </span><span>pos</span><span> </span><span>=</span><span> </span><span>4</span><span>)</span><span>
</span><span>}</span><span>
</span></code></pre>
</div>
<h4>Code to produce the histogram of the residuals</h4>
<div>
<pre><code><span>## Predicted points
</span><span>yhat</span><span> </span><span>&lt;-</span><span> </span><span>meany</span><span> </span><span>-</span><span> </span><span>r</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>meanx</span><span> </span><span>+</span><span> </span><span>r</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>x</span><span>
</span><span>## Prediction errors (residuals): actual - predicted
</span><span>ehat</span><span> </span><span>&lt;-</span><span> </span><span>y</span><span> </span><span>-</span><span> </span><span>yhat</span><span>
</span><span>hist</span><span>(</span><span>ehat</span><span>,</span><span> </span><span>breaks</span><span> </span><span>=</span><span> </span><span>&quot;Scott&quot;</span><span>,</span><span> </span><span>probability</span><span> </span><span>=</span><span> </span><span>TRUE</span><span>)</span><span>
</span></code></pre>
</div>
<h4>Code to produce the RMSE bands</h4>
<div>
<pre><code><span>## RMS error for regression
</span><span>
</span><span>plot</span><span>(</span><span>x</span><span>,</span><span> </span><span>y</span><span>,</span><span>
</span><span>xlim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>),</span><span> </span><span>ylim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>),</span><span> </span><span>col</span><span> </span><span>=</span><span> </span><span>&quot;lightgrey&quot;</span><span>,</span><span>
</span><span>xaxt</span><span> </span><span>=</span><span> </span><span>&quot;n&quot;</span><span>,</span><span> </span><span>yaxt</span><span> </span><span>=</span><span> </span><span>&quot;n&quot;</span><span>,</span><span> </span><span>xaxs</span><span> </span><span>=</span><span> </span><span>&quot;i&quot;</span><span>,</span><span> </span><span>yaxs</span><span> </span><span>=</span><span> </span><span>&quot;i&quot;</span><span>,</span><span>
</span><span>main</span><span> </span><span>=</span><span> </span><span>&quot;Regression line and 1 and 2 RMSE bands&quot;</span><span>,</span><span>
</span><span>xlab</span><span> </span><span>=</span><span> </span><span>&quot;Father height in inches&quot;</span><span>,</span><span> </span><span>ylab</span><span> </span><span>=</span><span> </span><span>&quot;Son height in inches&quot;</span><span>)</span><span>
</span><span>axp</span><span> </span><span>&lt;-</span><span> </span><span>seq</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>,</span><span> </span><span>by</span><span> </span><span>=</span><span> </span><span>2</span><span>)</span><span>
</span><span>axis</span><span>(</span><span>1</span><span>,</span><span> </span><span>at</span><span> </span><span>=</span><span> </span><span>axp</span><span>,</span><span> </span><span>labels</span><span> </span><span>=</span><span> </span><span>axp</span><span>)</span><span>
</span><span>axis</span><span>(</span><span>2</span><span>,</span><span> </span><span>at</span><span> </span><span>=</span><span> </span><span>axp</span><span>,</span><span> </span><span>labels</span><span> </span><span>=</span><span> </span><span>axp</span><span>)</span><span>
</span><span>## Regression line
</span><span>a</span><span> </span><span>&lt;-</span><span> </span><span>meany</span><span> </span><span>-</span><span> </span><span>r</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>meanx</span><span>
</span><span>b</span><span> </span><span>&lt;-</span><span> </span><span>r</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>a</span><span> </span><span>=</span><span> </span><span>a</span><span>,</span><span> </span><span>b</span><span> </span><span>=</span><span> </span><span>b</span><span>,</span><span> </span><span>lwd</span><span> </span><span>=</span><span> </span><span>4</span><span>,</span><span> </span><span>col</span><span> </span><span>=</span><span> </span><span>&quot;red&quot;</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>a</span><span> </span><span>=</span><span> </span><span>a</span><span> </span><span>-</span><span> </span><span>2</span><span> </span><span>*</span><span> </span><span>RMSE</span><span>,</span><span> </span><span>b</span><span> </span><span>=</span><span> </span><span>b</span><span>,</span><span> </span><span>lwd</span><span> </span><span>=</span><span> </span><span>1</span><span>,</span><span> </span><span>col</span><span> </span><span>=</span><span> </span><span>&quot;red&quot;</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>a</span><span> </span><span>=</span><span> </span><span>a</span><span> </span><span>-</span><span> </span><span>1</span><span> </span><span>*</span><span> </span><span>RMSE</span><span>,</span><span> </span><span>b</span><span> </span><span>=</span><span> </span><span>b</span><span>,</span><span> </span><span>lwd</span><span> </span><span>=</span><span> </span><span>1</span><span>,</span><span> </span><span>col</span><span> </span><span>=</span><span> </span><span>&quot;red&quot;</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>a</span><span> </span><span>=</span><span> </span><span>a</span><span> </span><span>+</span><span> </span><span>1</span><span> </span><span>*</span><span> </span><span>RMSE</span><span>,</span><span> </span><span>b</span><span> </span><span>=</span><span> </span><span>b</span><span>,</span><span> </span><span>lwd</span><span> </span><span>=</span><span> </span><span>1</span><span>,</span><span> </span><span>col</span><span> </span><span>=</span><span> </span><span>&quot;red&quot;</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>a</span><span> </span><span>=</span><span> </span><span>a</span><span> </span><span>+</span><span> </span><span>2</span><span> </span><span>*</span><span> </span><span>RMSE</span><span>,</span><span> </span><span>b</span><span> </span><span>=</span><span> </span><span>b</span><span>,</span><span> </span><span>lwd</span><span> </span><span>=</span><span> </span><span>1</span><span>,</span><span> </span><span>col</span><span> </span><span>=</span><span> </span><span>&quot;red&quot;</span><span>)</span><span>
</span></code></pre>
</div>
<h4>Code to produce the Graph of averages with strips of points</h4>
<div>
<pre><code><span>plot</span><span>(</span><span>x</span><span>,</span><span> </span><span>y</span><span>,</span><span>
</span><span>xlim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>),</span><span> </span><span>ylim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>),</span><span> </span><span>col</span><span>=</span><span>&quot;lightgrey&quot;</span><span>,</span><span>
</span><span>xaxt</span><span> </span><span>=</span><span> </span><span>&quot;n&quot;</span><span>,</span><span> </span><span>yaxt</span><span> </span><span>=</span><span> </span><span>&quot;n&quot;</span><span>,</span><span> </span><span>xaxs</span><span> </span><span>=</span><span> </span><span>&quot;i&quot;</span><span>,</span><span> </span><span>yaxs</span><span> </span><span>=</span><span> </span><span>&quot;i&quot;</span><span>,</span><span>
</span><span>main</span><span> </span><span>=</span><span> </span><span>&quot;Pearson's data. GoA and regression line&quot;</span><span>,</span><span>
</span><span>xlab</span><span> </span><span>=</span><span> </span><span>&quot;Father height in inches&quot;</span><span>,</span><span> </span><span>ylab</span><span> </span><span>=</span><span> </span><span>&quot;Son height in inches&quot;</span><span>)</span><span>
</span><span>axp</span><span> </span><span>&lt;-</span><span> </span><span>seq</span><span>(</span><span>58</span><span>,</span><span> </span><span>80</span><span>,</span><span> </span><span>by</span><span> </span><span>=</span><span> </span><span>2</span><span>)</span><span>
</span><span>axis</span><span>(</span><span>1</span><span>,</span><span> </span><span>at</span><span> </span><span>=</span><span> </span><span>axp</span><span>,</span><span> </span><span>labels</span><span> </span><span>=</span><span> </span><span>axp</span><span>)</span><span>
</span><span>axis</span><span>(</span><span>2</span><span>,</span><span> </span><span>at</span><span> </span><span>=</span><span> </span><span>axp</span><span>,</span><span> </span><span>labels</span><span> </span><span>=</span><span> </span><span>axp</span><span>)</span><span>
</span><span>## Point of averages (center of the cloud)
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>meanx</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;green&quot;</span><span>)</span><span>
</span><span>abline</span><span>(</span><span>h</span><span>=</span><span>meany</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;green&quot;</span><span>)</span><span>
</span><span>fround</span><span> </span><span>&lt;-</span><span> </span><span>66</span><span>
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>fround</span><span>+</span><span>c</span><span>(</span><span>-.5</span><span>,</span><span>.5</span><span>),</span><span> </span><span>lty</span><span>=</span><span>3</span><span>)</span><span>
</span><span>with</span><span>(</span><span>subset</span><span>(</span><span>father.son</span><span>,</span><span> </span><span>round</span><span>(</span><span>fheight</span><span>,</span><span>0</span><span>)</span><span> </span><span>==</span><span> </span><span>fround</span><span>),</span><span>
</span><span>points</span><span>(</span><span>fheight</span><span>,</span><span> </span><span>sheight</span><span>,</span><span> </span><span>pch</span><span>=</span><span>16</span><span>,</span><span> </span><span>col</span><span>=</span><span>ifelse</span><span>(</span><span>sheight</span><span> </span><span>></span><span> </span><span>meany</span><span> </span><span>+</span><span> </span><span>(</span><span>fheight</span><span>-</span><span>meanx</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>*</span><span>r</span><span>,</span><span> </span><span>&quot;darkgreen&quot;</span><span>,</span><span> </span><span>&quot;blue&quot;</span><span>)))</span><span>
</span><span>fround</span><span> </span><span>&lt;-</span><span> </span><span>70</span><span>
</span><span>abline</span><span>(</span><span>v</span><span>=</span><span>fround</span><span>+</span><span>c</span><span>(</span><span>-.5</span><span>,</span><span>.5</span><span>),</span><span> </span><span>lty</span><span>=</span><span>3</span><span>)</span><span>
</span><span>with</span><span>(</span><span>subset</span><span>(</span><span>father.son</span><span>,</span><span> </span><span>round</span><span>(</span><span>fheight</span><span>,</span><span>0</span><span>)</span><span> </span><span>==</span><span> </span><span>fround</span><span>),</span><span>
</span><span>points</span><span>(</span><span>fheight</span><span>,</span><span> </span><span>sheight</span><span>,</span><span>
</span><span>pch</span><span>=</span><span>16</span><span>,</span><span> </span><span>col</span><span>=</span><span>ifelse</span><span>(</span><span>sheight</span><span> </span><span>></span><span> </span><span>meany</span><span> </span><span>+</span><span> </span><span>(</span><span>fheight</span><span>-</span><span>meanx</span><span>)</span><span>/</span><span>sd</span><span>(</span><span>x</span><span>)</span><span>*</span><span>sd</span><span>(</span><span>y</span><span>)</span><span>*</span><span>r</span><span>,</span><span> </span><span>&quot;darkgreen&quot;</span><span>,</span><span> </span><span>&quot;blue&quot;</span><span>)))</span><span>
</span><span>## Graph of averages.
</span><span>sgav</span><span> </span><span>&lt;-</span><span> </span><span>with</span><span>(</span><span>father.son</span><span>,</span><span> </span><span>tapply</span><span>(</span><span>sheight</span><span>,</span><span> </span><span>round</span><span>(</span><span>fheight</span><span>,</span><span>0</span><span>),</span><span> </span><span>mean</span><span>))</span><span>
</span><span>## sgavnum &lt;- with(father.son, tapply(sheight, round(fheight,0), length))
</span><span>points</span><span>(</span><span>as.numeric</span><span>(</span><span>names</span><span>(</span><span>sgav</span><span>)),</span><span> </span><span>sgav</span><span>,</span><span> </span><span>col</span><span>=</span><span>&quot;red&quot;</span><span>,</span><span> </span><span>pch</span><span>=</span><span>16</span><span>)</span><span>
</span><span>## text(as.numeric(names(sgav)), sgav, sgavnum, pos=3)
</span></code></pre>
</div>
<h4>Code to produce the histograms of the strips</h4>
<div>
<pre><code><span>par</span><span>(</span><span>mfrow</span><span>=</span><span>c</span><span>(</span><span>1</span><span>,</span><span>2</span><span>))</span><span>
</span><span>fround</span><span> </span><span>&lt;-</span><span> </span><span>66</span><span>
</span><span>ystrip66</span><span> </span><span>&lt;-</span><span> </span><span>with</span><span>(</span><span>subset</span><span>(</span><span>father.son</span><span>,</span><span> </span><span>round</span><span>(</span><span>fheight</span><span>,</span><span>0</span><span>)</span><span> </span><span>==</span><span> </span><span>fround</span><span>),</span><span>
</span><span>sheight</span><span>)</span><span>
</span><span>hist</span><span>(</span><span>ystrip66</span><span>,</span><span> </span><span>breaks</span><span> </span><span>=</span><span> </span><span>7</span><span>)</span><span>
</span><span>fround</span><span> </span><span>&lt;-</span><span> </span><span>70</span><span>
</span><span>ystrip70</span><span> </span><span>&lt;-</span><span> </span><span>with</span><span>(</span><span>subset</span><span>(</span><span>father.son</span><span>,</span><span> </span><span>round</span><span>(</span><span>fheight</span><span>,</span><span>0</span><span>)</span><span> </span><span>==</span><span> </span><span>fround</span><span>),</span><span>
</span><span>sheight</span><span>)</span><span>
</span><span>hist</span><span>(</span><span>ystrip70</span><span>,</span><span> </span><span>breaks</span><span> </span><span>=</span><span> </span><span>7</span><span>)</span><span>
</span></code></pre>
</div>
<h3>Previous</h3>
<p><a target="_blank" rel="nofollow" href="https://rbertolusso.github.io/posts/LR01-correlation">LR01: Correlation</a></p>
<p><a target="_blank" rel="nofollow" href="https://rbertolusso.github.io/posts/LR02-SDline-GoA-regression">LR02: SD line, GoA, Regression</a></p>
<h3>Related</h3>
<p><a target="_blank" rel="nofollow" href="https://rbertolusso.github.io/posts/intubate-and-stat-functions-in-pipelines">intubate &lt;||> R stat functions in data science pipelines</a></p>
<p><div><div><div><div data-width="55" data-show-faces="false" data-send="true" data-layout="box_count" data-height="62" data-href="https://www.r-bloggers.com/lr03-residuals-and-rmse/"></div></div><div><a data-url="https://www.r-bloggers.com/lr03-residuals-and-rmse/" data-counturl="https://www.r-bloggers.com/lr03-residuals-and-rmse/" data-count="vertical" data-text="LR03: Residuals and RMSE" class="twitter-share-button" data-via="rbloggers" href="https://twitter.com/share"></a></div><div></div></div><div></div></div>
<div>
<div>To <strong>leave a comment</strong> for the author, please follow the link and comment on their blog: <strong><a href="http://rbertolusso.github.io/posts/LR03-residuals-RMSE"> intubate &lt;||>, XBRL, bioPN, sbioPN, and stats with R</a></strong>.</div>
<hr>
<a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers.com</a> offers <strong><a rel="nofollow" href="https://feedburner.google.com/fb/a/mailverify?uri=RBloggers">daily e-mail updates</a></strong> about <a title="The R Project for Statistical Computing" rel="nofollow" href="https://www.r-project.org/">R</a> news and <a title="R tutorials" rel="nofollow" href="https://www.r-bloggers.com/?s=tutorial">tutorials</a> on topics such as: <a title="Data science" rel="nofollow" href="https://www.r-bloggers.com/?s=data%20science">Data science</a>, <a title="Big Data" rel="nofollow" href="https://www.r-bloggers.com/?s=Big%20Data">Big Data, <a title="R jobs" rel="nofollow" href="https://www.r-users.com/">R jobs</a>, visualization (<a title="ggplot and ggplot2 tutorials" rel="nofollow" href="https://www.r-bloggers.com/?s=ggplot2">ggplot2</a>, <a title="Boxplots using lattice and ggplot2 tutorials" rel="nofollow" href="https://www.r-bloggers.com/?s=boxplot">Boxplots</a>, <a title="Maps and gis" rel="nofollow" href="https://www.r-bloggers.com/?s=map">maps</a>, <a title="Animation in R" rel="nofollow" href="https://www.r-bloggers.com/?s=animation">animation</a>), programming (<a title="RStudio IDE for R" rel="nofollow" href="https://www.r-bloggers.com/?s=RStudio">RStudio</a>, <a title="Sweave and literate programming" rel="nofollow" href="https://www.r-bloggers.com/?s=sweave">Sweave</a>, <a title="LaTeX in R" rel="nofollow" href="https://www.r-bloggers.com/?s=LaTeX">LaTeX</a>, <a title="SQL and databases" rel="nofollow" href="https://www.r-bloggers.com/?s=SQL">SQL</a>, <a title="Eclipse IDE for R" rel="nofollow" href="https://www.r-bloggers.com/?s=eclipse">Eclipse</a>, <a title="git and github, Version Control System" rel="nofollow" href="https://www.r-bloggers.com/?s=git">git</a>, <a title="Large data in R using Hadoop" rel="nofollow" href="https://www.r-bloggers.com/?s=hadoop">hadoop</a>, <a title="Web Scraping of google, facebook, yahoo, twitter and more using R" rel="nofollow" href="https://www.r-bloggers.com/?s=Web+Scraping">Web Scraping</a>) statistics (<a title="Regressions and ANOVA analysis tutorials" rel="nofollow" href="https://www.r-bloggers.com/?s=regression">regression</a>, <a title="principal component analysis tutorial" rel="nofollow" href="https://www.r-bloggers.com/?s=PCA">PCA</a>, <a title="Time series" rel="nofollow" href="https://www.r-bloggers.com/?s=time+series">time series</a>, <a title="finance trading" rel="nofollow" href="https://www.r-bloggers.com/?s=trading">trading</a>) and more...
<img height="1" alt="" width="1" src="http://feeds.feedburner.com/~r/RBloggers/~4/fwVcTiVFI1s"></a></div>