<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/principal-component-analysis/"></div></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="http://www.aaronschlegel.com/principal-component-analysis/"> R – Aaron Schlegel</a></strong>, and kindly contributed to <a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers)</a> </div>
<p>Often, it is not helpful or informative to only look at all the variables in a dataset for correlations or covariances. A preferable approach is to derive new variables from the original variables that preserve most of the information given by their variances. Principal component analysis is a widely used and popular statistical method for reducing data with many dimensions (variables) by projecting the data with fewer dimensions using linear combinations of the variables, known as principal components. The new projected variables (principal components) are uncorrelated with each other and are ordered so that the first few components retain most of the variation present in the original variables. Thus, PCA is also useful in situations where the independent variables are correlated with each other and can be employed in exploratory data analysis or for making predictive models. Principal component analysis can also reveal important features of the data such as outliers and departures from a multinormal distribution.</p>
<h6>Defining Principal Components</h6>
<p>The first step in defining the principal components of <span data-display="false" class="wp-katex-eq">p</span> original variables is to find a linear function <span data-display="false" class="wp-katex-eq">a_1’y</span>, where <span data-display="false" class="wp-katex-eq">a_1</span> is a vector of <span data-display="false" class="wp-katex-eq">p</span> constants, for the observation vectors that have maximum variance. This linear function is defined as:</p>
<p><span data-display="true" class="wp-katex-eq"> a_1’y = a_{11}x_1 + a_{12}x_2 + \cdots + a_{1p}x_p = \sum_{j=1}^p a_{1j}x_j </span></p>
<p>Principal component analysis continues to find a linear function <span data-display="false" class="wp-katex-eq">a_2’y</span> that is uncorrelated with <span data-display="false" class="wp-katex-eq">a_1’y</span> with maximized variance and so on up to <span data-display="false" class="wp-katex-eq">k</span> principal components.</p>
<h6>Derivation of Principal Components</h6>
<p>The principal components of a dataset are obtained from the sample covariance matrix <span data-display="false" class="wp-katex-eq">S</span> or the correlation matrix <span data-display="false" class="wp-katex-eq">R</span>. Although principal components obtained from <span data-display="false" class="wp-katex-eq">S</span> is the original method of principal component analysis, components from <span data-display="false" class="wp-katex-eq">R</span> may be more interpretable if the original variables have different units or wide variances (Rencher 2002, pp. 393). For now, <span data-display="false" class="wp-katex-eq">S</span> will be referred to as <span data-display="false" class="wp-katex-eq">\Sigma</span> (denotes a known covariance matrix) which will be used in the derivation.</p>
<p>The goal of the derivation is to find <span data-display="false" class="wp-katex-eq">a’_ky</span> that maximizes the variance of <span data-display="false" class="wp-katex-eq">a’_ky \Sigma a_k</span>. For this, we will consider the first vector <span data-display="false" class="wp-katex-eq">a’_1y</span> that maximizes <span data-display="false" class="wp-katex-eq">Var(a’_1y) = a’_1y \Sigma a_1</span>. To do this maximization, we will need a constraint to rein in otherwise unnecessarily large values of <span data-display="false" class="wp-katex-eq">a_1</span>. The constraint in this example is the unit length vector <span data-display="false" class="wp-katex-eq">a_1′ a_1 = 1</span>. This constraint is employed with a Lagrange multiplier <span data-display="false" class="wp-katex-eq">\lambda</span> so that the function is maximized at an equality constraint of <span data-display="false" class="wp-katex-eq">g(x) = 0</span>. Thus the Lagrangian function is defined as:</p>
<p><span data-display="true" class="wp-katex-eq"> a’_1 \Sigma a_1 – \lambda(a_1’y a_1 – 1) </span></p>
<h6>Brief Aside: Lagrange Multipliers</h6>
<p>The Lagrange mulitiplier method is used for finding a maximum or minimum of a multivariate function with some constraint on the input values. As we are interested in maximization, the problem can be briefly stated as ’maximize <span data-display="false" class="wp-katex-eq">f(x)</span> subject to <span data-display="false" class="wp-katex-eq">g(x) = c</span>’. In this example, <span data-display="false" class="wp-katex-eq">g(x) = a_1’y a_1 = 1</span> and <span data-display="false" class="wp-katex-eq">f(x) = a_1’y \Sigma a_1</span>. The Lagrange multiplier, defined as <span data-display="false" class="wp-katex-eq">\lambda</span> allows the combination of <span data-display="false" class="wp-katex-eq">f(x)</span> and <span data-display="false" class="wp-katex-eq">g(x)</span> into a new function <span data-display="false" class="wp-katex-eq">L(x, \lambda)</span>, defined as:</p>
<p><span data-display="true" class="wp-katex-eq"> L(a_1’y, \lambda) = f(a_1’y) – \lambda(g(a_1’y) – c) </span></p>
<p>The sign of <span data-display="false" class="wp-katex-eq">\lambda</span> can be positive or negative. The new function is then solved for a stationary point, in this case <span data-display="false" class="wp-katex-eq">0</span>, using partial derivatives:</p>
<p><span data-display="true" class="wp-katex-eq"> \frac{\partial L(a_1’y, \lambda)}{\partial a_1’y} = 0 </span></p>
<p>Returning to principal component analysis, we differentiate <span data-display="false" class="wp-katex-eq">L(a_1) = a’_1 \Sigma a_1 – \lambda(a_1’y a_1 – 1)</span> with respect to <span data-display="false" class="wp-katex-eq">a_1</span>:</p>
<p><span data-display="true" class="wp-katex-eq"> \frac{\partial L}{\partial a_1} = 2\Sigma a_1 – 2\lambda a_1 = 0 </span><br>
<span data-display="true" class="wp-katex-eq"> \Sigma a_1 – \lambda a_1 = 0 </span></p>
<p>Expressing the above with an identity matrix, <span data-display="false" class="wp-katex-eq">I</span>:</p>
<p><span data-display="true" class="wp-katex-eq"> (\Sigma – \lambda I) a_1 = 0 </span></p>
<p>Which shows <span data-display="false" class="wp-katex-eq">\lambda</span> is an eigenvector of the covariance matrix <span data-display="false" class="wp-katex-eq">\Sigma</span> and <span data-display="false" class="wp-katex-eq">a_1</span> is the corresponding eigenvector. As stated previously, we are interested in finding <span data-display="false" class="wp-katex-eq">a_1’y</span> with maximum variance. Therefore <span data-display="false" class="wp-katex-eq">\lambda</span> must be as large as possible which follows <span data-display="false" class="wp-katex-eq">a_1</span> is the eigenvector corresponding to the largest eigenvalue of <span data-display="false" class="wp-katex-eq">\Sigma</span>.</p>
<p>The remaining principal components are found in a similar manner and correspond to the <span data-display="false" class="wp-katex-eq">k</span>th principal component. Thus the second principal component is <span data-display="false" class="wp-katex-eq">a_2’y</span> and is equivalent to the eigenvector of the second largest eigenvalue of <span data-display="false" class="wp-katex-eq">\Sigma</span>, and so on.</p>
<h6>Principal Component Analysis</h6>
<p>Twenty engineer apprentices and twenty pilots were given six tests. The data were obtained from the <a target="_blank" rel="nofollow">companion FTP site</a> of the book Methods of Multivariate Analysis by Alvin Rencher. The tests measured the following attributes:</p>
<ul>
<li>Intelligence</li>
<li>Form Relations</li>
<li>Dynamometer</li>
<li>Dotting</li>
<li>Sensory Motor Coordination</li>
<li>Perservation</li>
</ul>
<p>Principal component analysis will be performed on the data to transform the attributes into new variables that will hopefully be more open to interpretation and allow us to find any irregularities in the data such as outliers.</p>
<p>Load the data and name the columns. The factors in the <code>Group</code> column are renamed to their actual grouping names.</p>
<pre>
pilots &lt;- read.table('PILOTS.DAT', col.names =  c('Group', 'Intelligence', 'Form Relations',
                                                  'Dynamometer', 'Dotting', 'Sensory Motor Coordination',
                                                  'Perservation'))
pilots$Group &lt;- ifelse(pilots$Group == 1, 'Apprentice', 'Pilot')
</pre>
<p>Inspect the first few rows of the data.</p>
<pre>
head(pilots)
##        Group Intelligence Form.Relations Dynamometer Dotting
## 1 Apprentice          121             22          74     223
## 2 Apprentice          108             30          80     175
## 3 Apprentice          122             49          87     266
## 4 Apprentice           77             37          66     178
## 5 Apprentice          140             35          71     175
## 6 Apprentice          108             37          57     241
##   Sensory.Motor.Coordination Perservation
## 1                         54          254
## 2                         40          300
## 3                         41          223
## 4                         80          209
## 5                         38          261
## 6                         59          245
</pre>
<p>The variables appear to be measured in different units which may lead to the variables with larger variances dominating the principal components of the covariance matrix <span data-display="false" class="wp-katex-eq">S</span>. We will perform principal component analysis on the correlation matrix <span data-display="false" class="wp-katex-eq">R</span> later in the example to find a scaled and more balanced representation of the components.</p>
<p>Find the covariance matrix <span data-display="false" class="wp-katex-eq">S</span> of the data. The grouping column is not included.</p>
<pre>
S &lt;- cov(pilots[,2:7])
S
##                            Intelligence Form.Relations Dynamometer
## Intelligence                  528.19487       35.98974    27.97949
## Form.Relations                 35.98974       68.91282   -12.09744
## Dynamometer                    27.97949      -12.09744   145.18974
## Dotting                       104.42821      -81.75128   128.88205
## Sensory.Motor.Coordination    -20.03077      -33.00513   -30.85641
## Perservation                  291.15385      -18.28205    29.38462
##                               Dotting Sensory.Motor.Coordination
## Intelligence                104.42821                  -20.03077
## Form.Relations              -81.75128                  -33.00513
## Dynamometer                 128.88205                  -30.85641
## Dotting                    1366.43013                 -113.58077
## Sensory.Motor.Coordination -113.58077                  264.35641
## Perservation                395.18590                  -79.85897
##                            Perservation
## Intelligence                  291.15385
## Form.Relations                -18.28205
## Dynamometer                    29.38462
## Dotting                       395.18590
## Sensory.Motor.Coordination    -79.85897
## Perservation                 1069.11538
</pre>
<p>The total variance is defined as:</p>
<p><span data-display="true" class="wp-katex-eq"> \sum^k_{j=1} s_{jj} </span></p>
<p>Which is also equal to the sum of the eigenvalues of <span data-display="false" class="wp-katex-eq">S</span>.</p>
<pre>
sum(diag(S))
## [1] 3442.199
</pre>
<p>Compute the eigenvalues and corresponding eigenvectors of <span data-display="false" class="wp-katex-eq">S</span>.</p>
<pre>
s.eigen &lt;- eigen(S)
s.eigen
## $values
## [1] 1722.0424  878.3578  401.4386  261.0769  128.9051   50.3785
## 
## $vectors
##             [,1]        [,2]        [,3]        [,4]        [,5]
## [1,] -0.21165160 -0.38949336  0.88819049  0.03082062 -0.04760343
## [2,]  0.03883125 -0.06379320  0.09571590 -0.19128493 -0.14793191
## [3,] -0.08012946  0.06602004  0.08145863 -0.12854488  0.97505667
## [4,] -0.77552673  0.60795970  0.08071120  0.08125631 -0.10891968
## [5,]  0.09593926 -0.01046493  0.01494473  0.96813856  0.10919120
## [6,] -0.58019734 -0.68566916 -0.43426141  0.04518327  0.03644629
##             [,6]
## [1,]  0.10677164
## [2,] -0.96269790
## [3,] -0.12379748
## [4,] -0.06295166
## [5,] -0.20309559
## [6,] -0.03572141
</pre>
<p>The eigenvectors represent the principal components of <span data-display="false" class="wp-katex-eq">S</span>. The eigenvalues of <span data-display="false" class="wp-katex-eq">S</span> are used to find the proportion of the total variance explained by the components.</p>
<pre>
for (s in s.eigen$values) {
  print(s / sum(s.eigen$values))
}
## [1] 0.5002739
## [1] 0.2551734
## [1] 0.1166227
## [1] 0.07584597
## [1] 0.03744848
## [1] 0.01463556
</pre>
<p>The first two principal components account for 75.5% of the total variance. A scree graph of the eigenvalues can be plotted to visualize the proportion of variance explained by each subsequential eigenvalue.</p>
<pre>
plot(s.eigen$values, xlab = 'Eigenvalue Number', ylab = 'Eigenvalue Size', main = 'Scree Graph')
lines(s.eigen$values)
</pre>
<p><a target="_blank" rel="nofollow" href="http://www.aaronschlegel.com/?attachment_id=22448"><img sizes="(max-width: 1060px) 100vw, 1060px" alt="Scree plot from principal component analysis" class="wp-image-22448" src="https://i0.wp.com/www.aaronschlegel.com/wp-content/uploads/scree_plot.png?resize=1060%2C757" data-recalc-dims="1"></a></p>
<p>The elements of the eigenvectors of <span data-display="false" class="wp-katex-eq">S</span> are the ‘coefficients’ or ‘loadings’ of the principal components.</p>
<pre>
s.eigen$vectors
##             [,1]        [,2]        [,3]        [,4]        [,5]
## [1,] -0.21165160 -0.38949336  0.88819049  0.03082062 -0.04760343
## [2,]  0.03883125 -0.06379320  0.09571590 -0.19128493 -0.14793191
## [3,] -0.08012946  0.06602004  0.08145863 -0.12854488  0.97505667
## [4,] -0.77552673  0.60795970  0.08071120  0.08125631 -0.10891968
## [5,]  0.09593926 -0.01046493  0.01494473  0.96813856  0.10919120
## [6,] -0.58019734 -0.68566916 -0.43426141  0.04518327  0.03644629
##             [,6]
## [1,]  0.10677164
## [2,] -0.96269790
## [3,] -0.12379748
## [4,] -0.06295166
## [5,] -0.20309559
## [6,] -0.03572141
</pre>
<p>The first two principal components are thus:</p>
<p><span data-display="true" class="wp-katex-eq"> z_1 = a’_1y = -.212y_1 + .039y_2 – 0.080y_3 – 0.776y_4 + 0.096y_5 – 0.580y_6 </span> <span data-display="true" class="wp-katex-eq"> z_2 = a’_2y = -.389y_1 – .064y_2 + 0.066y_3 + 0.608y_4 – 0.010y_5 – 0.686y_6 </span></p>
<h6>Principal Component Analysis with R</h6>
<p>Computing the principal components in R is straightforward with the functions <code>prcomp()</code> and <code>princomp()</code>. The difference between the two is simply the method employed to calculate PCA. According to <code>?prcomp</code>:</p>
<blockquote><p>
The calculation is done by a singular value decomposition of the (centered and possibly scaled) data matrix, not by using eigen on the covariance matrix. This is generally the preferred method for numerical accuracy.
</p>
</blockquote>
<p>From <code>?princomp</code>:</p>
<blockquote><p>
The calculation is done using eigen on the correlation or covariance matrix, as determined by cor. This is done for compatibility with the S-PLUS result.
</p>
</blockquote>
<pre>
pilots.pca &lt;- prcomp(pilots[,2:7])
pilots.pca
## Standard deviations:
## [1] 41.497499 29.637102 20.035932 16.157875 11.353640  7.097781
## 
## Rotation:
##                                    PC1         PC2         PC3         PC4
## Intelligence                0.21165160 -0.38949336  0.88819049 -0.03082062
## Form.Relations             -0.03883125 -0.06379320  0.09571590  0.19128493
## Dynamometer                 0.08012946  0.06602004  0.08145863  0.12854488
## Dotting                     0.77552673  0.60795970  0.08071120 -0.08125631
## Sensory.Motor.Coordination -0.09593926 -0.01046493  0.01494473 -0.96813856
## Perservation                0.58019734 -0.68566916 -0.43426141 -0.04518327
##                                    PC5         PC6
## Intelligence               -0.04760343 -0.10677164
## Form.Relations             -0.14793191  0.96269790
## Dynamometer                 0.97505667  0.12379748
## Dotting                    -0.10891968  0.06295166
## Sensory.Motor.Coordination  0.10919120  0.20309559
## Perservation                0.03644629  0.03572141
</pre>
<p>Although we didn’t use the preferred method of applying singular value decomposition, the components reported by the <code>prcomp()</code> are the same as what was computed earlier save arbitrary scalings of <span data-display="false" class="wp-katex-eq">-1</span> to some of the eigenvectors.</p>
<p>The summary method of <code>prcomp()</code> also outputs the proportion of variance explained by the components.</p>
<pre>
summary(pilots.pca)
## Importance of components:
##                            PC1     PC2     PC3      PC4      PC5     PC6
## Standard deviation     41.4975 29.6371 20.0359 16.15788 11.35364 7.09778
## Proportion of Variance  0.5003  0.2552  0.1166  0.07585  0.03745 0.01464
## Cumulative Proportion   0.5003  0.7554  0.8721  0.94792  0.98536 1.00000
</pre>
<h6>Plotting of Principal Components</h6>
<p>The first two principal components are often plotted as a scatterplot which may reveal interesting features of the data, such as departures from normality, outliers or non-linearity. The first two principal components are evaluated for each observation vector and plotted.</p>
<p>The <a target="_blank" rel="nofollow" href="https://cran.r-project.org/web/packages/ggfortify/index.html">ggfortify package</a> provides a handy method for plotting the first two principal components with <code>autoplot()</code>.</p>
<pre>
library(ggfortify)
</pre>
<p>The <code>autoplot()</code> function also generates a useful data table of the calculated principal components we which we will use later.</p>
<pre>
pca.plot &lt;- autoplot(pilots.pca, data = pilots, colour = 'Group')
pca.plot
</pre>
<p><a target="_blank" rel="nofollow" href="http://www.aaronschlegel.com/?attachment_id=22445"><img sizes="(max-width: 1060px) 100vw, 1060px" alt="plot of first two principal components" class="wp-image-22445" src="https://i2.wp.com/www.aaronschlegel.com/wp-content/uploads/pc_plot1.png?resize=1060%2C757" data-recalc-dims="1"></a></p>
<p>The points of the two groups are clustered for the most part; however, the three points at the top of the graph may be outliers. The data does not appear to depart widely from multivariate normality. We will see if this conclusion changes when the PCs from the correlation matrix <span data-display="false" class="wp-katex-eq">R</span> are plotted.</p>
<p>To recreate the graph generated by <code>autoplot()</code>, scale the data using the standard deviations of the principal components multiplied by the square root of the number of observations. The principal components are then computed for each observation vector. Note the first eigenvector is multiplied by a scaling factor of <span data-display="false" class="wp-katex-eq">-1</span> so the signs what was reported by the <code>prcomp()</code> function.</p>
<pre>
scaling &lt;- pilots.pca$sdev[1:2] * sqrt(nrow(pilots))

pc1 &lt;- rowSums(t(t(sweep(pilots[,2:7], 2 ,colMeans(pilots[,2:7]))) * s.eigen$vectors[,1] * -1) / scaling[1])
pc2 &lt;- rowSums(t(t(sweep(pilots[,2:7], 2, colMeans(pilots[,2:7]))) * s.eigen$vectors[,2]) / scaling[2])
</pre>
<p>Collect the PCs in a <code>data.frame</code> and plot using <code>ggplot</code> (loaded when <code>ggfortify</code> was loaded).</p>
<pre>
df &lt;- data.frame(pc1, pc2, c(rep('Apprentice', 20), rep('Pilot', 20)))
colnames(df) &lt;- c('PC1', 'PC2', 'Group')

ggplot(df, aes(x=PC1, y=PC2, color=Group)) + 
  geom_point()
</pre>
<p><a target="_blank" rel="nofollow" href="http://www.aaronschlegel.com/?attachment_id=22446"><img sizes="(max-width: 1060px) 100vw, 1060px" alt="recreated plot of first two principal components" class="wp-image-22446" src="https://i2.wp.com/www.aaronschlegel.com/wp-content/uploads/pc_plot2.png?resize=1060%2C757" data-recalc-dims="1"></a></p>
<p>The scaling employed when calculating the PCs can be omitted. To remove scaling in the <code>autplot()</code> function, set the <code>scaling</code> argument to 0.</p>
<h6>Principal Component Analysis with the Correlation Matrix <span data-display="false" class="wp-katex-eq">R</span></h6>
<p>As mentioned previously, although principal component analysis is typically performed on the covariance matrix <span data-display="false" class="wp-katex-eq">S</span>, it often makes more intuitive sense to apply PCA to the correlation matrix. Cases where using <span data-display="false" class="wp-katex-eq">R</span> may be preferable to <span data-display="false" class="wp-katex-eq">S</span> include data that is measured in different units or has wide variances. The pilot data analyzed does not appear to have commensurate units for each variable, and because we have very little information regarding the tests and the measurements collected, it might make sense to employ the <span data-display="false" class="wp-katex-eq">R</span> matrix rather than <span data-display="false" class="wp-katex-eq">S</span>.</p>
<p>The correlation matrix is found with the <code>cor()</code> function.</p>
<pre>
R &lt;- cor(pilots[,2:7])
R
##                            Intelligence Form.Relations Dynamometer
## Intelligence                 1.00000000     0.18863907  0.10103566
## Form.Relations               0.18863907     1.00000000 -0.12094150
## Dynamometer                  0.10103566    -0.12094150  1.00000000
## Dotting                      0.12292126    -0.26641020  0.28935484
## Sensory.Motor.Coordination  -0.05360504    -0.24453244 -0.15750071
## Perservation                 0.38744776    -0.06735388  0.07458298
##                               Dotting Sensory.Motor.Coordination
## Intelligence                0.1229213                -0.05360504
## Form.Relations             -0.2664102                -0.24453244
## Dynamometer                 0.2893548                -0.15750071
## Dotting                     1.0000000                -0.18898014
## Sensory.Motor.Coordination -0.1889801                 1.00000000
## Perservation                0.3269606                -0.15021611
##                            Perservation
## Intelligence                 0.38744776
## Form.Relations              -0.06735388
## Dynamometer                  0.07458298
## Dotting                      0.32696061
## Sensory.Motor.Coordination  -0.15021611
## Perservation                 1.00000000
</pre>
<p>Find the eigenvalues and eigenvectors of the <span data-display="false" class="wp-katex-eq">R</span> matrix.</p>
<pre>
r.eigen &lt;- eigen(R)
</pre>
<p>As with the covariance matrix, we can compute the proportion of total variance explained by the eigenvalues.</p>
<pre>
for (r in r.eigen$values) {
  print(r / sum(r.eigen$values))
}
## [1] 0.2958546
## [1] 0.225736
## [1] 0.1787751
## [1] 0.1357993
## [1] 0.08843547
## [1] 0.07539955
</pre>
<p>What is readily noticeable is the first eigenvalue accounts for 30% of total variance compared with 50% of the variance of the <span data-display="false" class="wp-katex-eq">S</span> matrix. The first two components of <span data-display="false" class="wp-katex-eq">R</span> only account for 52% of total variance while the last two components have little significance. Thus, one may want to keep the first four components rather than the first two or three with the <span data-display="false" class="wp-katex-eq">S</span> matrix.</p>
<p>To perform principal component analysis using the correlation matrix using the <code>prcomp()</code> function, set the <code>scale</code> argument to <code>TRUE</code>.</p>
<pre>
pilots.pca.scaled &lt;- prcomp(pilots[,2:7], scale = TRUE)
pilots.pca.scaled
## Standard deviations:
## [1] 1.3323392 1.1637937 1.0356884 0.9026604 0.7284317 0.6726049
## 
## Rotation:
##                                    PC1        PC2        PC3        PC4
## Intelligence                0.40239072 -0.3964661  0.4617841 -0.3928149
## Form.Relations             -0.09715877 -0.7472294 -0.1752970 -0.1315611
## Dynamometer                 0.38541311  0.2181560 -0.4329575 -0.7177525
## Dotting                     0.54333623  0.3144601 -0.1065065  0.2453920
## Sensory.Motor.Coordination -0.31188931  0.3559400  0.6268314 -0.3992852
## Perservation                0.53629229 -0.1062657  0.4053555  0.3058981
##                                   PC5        PC6
## Intelligence               -0.2103062 -0.5187674
## Form.Relations             -0.2801896  0.5528697
## Dynamometer                 0.2585104  0.1855163
## Dotting                    -0.7066663  0.1869825
## Sensory.Motor.Coordination -0.2012981  0.4279773
## Perservation                0.5201339  0.4155385
</pre>
<pre>
summary(pilots.pca.scaled)
## Importance of components:
##                           PC1    PC2    PC3    PC4     PC5    PC6
## Standard deviation     1.3323 1.1638 1.0357 0.9027 0.72843 0.6726
## Proportion of Variance 0.2959 0.2257 0.1788 0.1358 0.08844 0.0754
## Cumulative Proportion  0.2959 0.5216 0.7004 0.8362 0.92460 1.0000
</pre>
<p>Plot the first two PCs of the correlation matrix using the <code>autoplot()</code> function.</p>
<pre>
pca.plot.scaled &lt;- autoplot(pilots.pca.scaled, data = pilots, colour = 'Group')
pca.plot.scaled
</pre>
<p><a target="_blank" rel="nofollow" href="http://www.aaronschlegel.com/?attachment_id=22447"><img sizes="(max-width: 1060px) 100vw, 1060px" alt="Plot of first two principal components obtained from correlation matrix" class="wp-image-22447" src="https://i2.wp.com/www.aaronschlegel.com/wp-content/uploads/pc_plot3.png?resize=1060%2C757" data-recalc-dims="1"></a></p>
<p>The points remain clustered for the most part; however, there does appear to be more points outside that may be considered outliers, though they don’t appear to be too far off from the cluster. It is important to remember the first two PCs of the <span data-display="false" class="wp-katex-eq">R</span> matrix only represent 52% of the total variance and thus may not be fully representative of the variance in the dataset.</p>
<h6>Interpreting Principal Components</h6>
<p>Interpretation of principal components is still a heavily researched topic in statistics, and although the components may be readily interpreted in most settings, this is not always the case (Joliffe, 2002).</p>
<p>One method of interpretation of the principal components is to calculate the correlation between the original data and the component. The <code>autoplot()</code> function also generates a nice data table with the original variables and the calculated PCs, which we will use here to find the correlations.</p>
<p>First, compute the correlations between the data and the calculated components of the covariance matrix <span data-display="false" class="wp-katex-eq">S</span>.</p>
<pre>
comps &lt;- pca.plot$data[,c(1:2,4:13)]
cor(comps[,3:8], comps[,c(1:2,9)])
##                                   PC1         PC2         PC3
## Intelligence                0.3821610 -0.50227170  0.77431660
## Form.Relations             -0.1941124 -0.22775091  0.23101683
## Dynamometer                 0.2759600  0.16238412  0.13544989
## Dotting                     0.8706127  0.48743512  0.04374714
## Sensory.Motor.Coordination -0.2448631 -0.01907555  0.01841632
## Perservation                0.7363518 -0.62149561 -0.26610222
</pre>
<p>The PCs can then be interpreted based on which variables they are most correlated in either a positive or negative direction. The level at which the correlations are significant is left to the researcher.</p>
<p>The first component is correlated negatively with Dotting, Perservation, Intelligence and Dynamometer. This correlation suggests the five variables vary together and when one goes down, the others decrease as well. The component is most correlated with Dotting at <span data-display="false" class="wp-katex-eq">-.087</span> and could be considered as primarily a measure of Dotting.</p>
<p>The second component is most correlated with Perservation and Intelligence, both in a negative direction. Dotting is correlated with the second component in a positive direction, which would indicate that as Perservation and Intelligence decrease, Dotting increases.</p>
<p>The third component is primarily correlated with Intelligence and not much else. This component could be viewed as a measure of the intelligence of the individual apprentice or pilot.</p>
<p>It was decided previously that due to lack of information regarding the variables and their units of measurement, it makes more sense to use the correlation matrix <span data-display="false" class="wp-katex-eq">R</span> for performing principal component analysis. Let’s see how the interpretation of the principal components changes when we use the <span data-display="false" class="wp-katex-eq">R</span> matrix.</p>
<pre>
comps.scaled &lt;- pca.plot.scaled$data[,c(1:2,4:13)]
cor(comps.scaled[,3:8], comps.scaled[,c(1:2,9)])
##                                   PC1        PC2        PC3
## Intelligence                0.5361209 -0.4614047  0.4782644
## Form.Relations             -0.1294484 -0.8696209 -0.1815530
## Dynamometer                 0.5135010  0.2538886 -0.4484090
## Dotting                     0.7239081  0.3659667 -0.1103075
## Sensory.Motor.Coordination -0.4155423  0.4142408  0.6492020
## Perservation                0.7145232 -0.1236713  0.4198220
</pre>
<p>The most apparent changes between the correlations of the original variables and the PCs of the <span data-display="false" class="wp-katex-eq">S</span> and <span data-display="false" class="wp-katex-eq">R</span> matrices are in components 2 and 3. The first principal component is still strongly correlated with the variables Dotting and Perservation, but now the variables Intelligence and Dynamometer are much more correlated and could indicate that as the former two variables decrease, the latter two increase.</p>
<p>The second component is now correlated the most with Forming Relations and not much else, whereas with the <span data-display="false" class="wp-katex-eq">S</span> matrix the component was correlated more to Perservation and Intelligence. This difference in variable correlations between the components of the two matrices may indicate Perservation and Intelligence were unduly dominating the variances.</p>
<p>The third component is now most correlated with Sensory Motor Coordination and secondarily Intelligence and Perservation, which indicates that subjects with high Sensory Motor Coordination test scores also have higher Intelligence and Perservation scores.</p>
<h6>Summary</h6>
<p>This post ended up being much longer than I had anticipated but I hope it is a good introduction to the power and benefits of principal component analysis. The post covered PCA with the covariance and correlation matrices as well as plotting and interpreting the principal components. I plan to continue discussing PCA in the future as there are many more topics and applications related to the dimension reduction technique.</p>
<h6>References</h6>
<p>Joliffe, I. T. (2002). Principal Component Analysis (2nd ed.). Springer.</p>
<p>Rencher, A. (n.d.). Methods of Multivariate Analysis (2nd ed.)</p>
<p><a target="_blank" rel="nofollow" href="https://onlinecourses.science.psu.edu/stat505/node/54">https://onlinecourses.science.psu.edu/stat505/node/54</a></p>
<p>The post <a target="_blank" rel="nofollow" href="http://www.aaronschlegel.com/principal-component-analysis/">Principal Component Analysis</a> appeared first on <a target="_blank" rel="nofollow" href="http://www.aaronschlegel.com/">Aaron Schlegel</a>.</p>
<p>
<div>
<div>To <strong>leave a comment</strong> for the author, please follow the link and comment on their blog: <strong><a href="http://www.aaronschlegel.com/principal-component-analysis/"> R – Aaron Schlegel</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/search/tutorial">tutorials</a> on topics such as: <a title="Data science" rel="nofollow" href="https://www.r-bloggers.com/search/data%20science">Data science</a>, <a title="Big Data" rel="nofollow" href="https://www.r-bloggers.com/search/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/search/ggplot2">ggplot2</a>, <a title="Boxplots using lattice and ggplot2 tutorials" rel="nofollow" href="https://www.r-bloggers.com/search/boxplot">Boxplots</a>, <a title="Maps and gis" rel="nofollow" href="https://www.r-bloggers.com/search/map">maps</a>, <a title="Animation in R" rel="nofollow" href="https://www.r-bloggers.com/search/animation">animation</a>), programming (<a title="RStudio IDE for R" rel="nofollow" href="https://www.r-bloggers.com/search/RStudio">RStudio</a>, <a title="Sweave and literate programming" rel="nofollow" href="https://www.r-bloggers.com/search/sweave">Sweave</a>, <a title="LaTeX in R" rel="nofollow" href="https://www.r-bloggers.com/search/LaTeX">LaTeX</a>, <a title="SQL and databases" rel="nofollow" href="https://www.r-bloggers.com/search/SQL">SQL</a>, <a title="Eclipse IDE for R" rel="nofollow" href="https://www.r-bloggers.com/search/eclipse">Eclipse</a>, <a title="git and github, Version Control System" rel="nofollow" href="https://www.r-bloggers.com/search/git">git</a>, <a title="Large data in R using Hadoop" rel="nofollow" href="https://www.r-bloggers.com/search/hadoop">hadoop</a>, <a title="Web Scraping of google, facebook, yahoo, twitter and more using R" rel="nofollow" href="https://www.r-bloggers.com/search/Web+Scraping">Web Scraping</a>) statistics (<a title="Regressions and ANOVA analysis tutorials" rel="nofollow" href="https://www.r-bloggers.com/search/regression">regression</a>, <a title="principal component analysis tutorial" rel="nofollow" href="https://www.r-bloggers.com/search/PCA">PCA</a>, <a title="Time series" rel="nofollow" href="https://www.r-bloggers.com/search/time+series">time series</a>, <a title="finance trading" rel="nofollow" href="https://www.r-bloggers.com/search/trading">trading</a>) and more...
<img height="1" alt="" width="1" src="http://feeds.feedburner.com/~r/RBloggers/~4/P2-97wmbo5k"></a></div>