# Modern Dimension Reduction

Philip D. Waggoner
University of Chicago
<EMAIL_ADDRESS>
###### Abstract

Data are not only ubiquitous in society, but are increasingly complex both in
size and dimensionality. Dimension reduction offers researchers and scholars
the ability to make such complex, high dimensional data spaces simpler and
more manageable. This Element offers readers a suite of modern unsupervised
dimension reduction techniques along with hundreds of lines of R code, to
efficiently represent the original high dimensional data space in a
simplified, lower dimensional subspace. Launching from the earliest dimension
reduction technique principal components analysis and using real social
science data, I introduce and walk readers through application of the
following techniques: locally linear embedding, t-distributed stochastic
neighbor embedding (t-SNE), uniform manifold approximation and projection,
self-organizing maps, and deep autoencoders. The result is a well-stocked
toolbox of unsupervised algorithms for tackling the complexities of high
dimensional data so common in modern society. All code is publicly accessible
on Github.

Forthcoming, Cambridge University Press000Please let me know if you plan to
cite this book in published work.

###### Contents

  1. 1 Introduction
    1. 1.1 Defining the Title
    2. 1.2 Running Example: 2019 American National Election Pilot Study
    3. 1.3 The Methods
    4. 1.4 The Methods Not Covered
    5. 1.5 Tidy Programming in R
    6. 1.6 Cleaning the Data
    7. 1.7 Addressing Missing Data
  2. 2 A Classic Approach to Dimension Reduction
    1. 2.1 Why PCA?
    2. 2.2 What is PCA Doing?
    3. 2.3 Formalizing PCA
    4. 2.4 Applying PCA to the ANES Data
    5. 2.5 Suggestions for Further Reading
  3. 3 Locally Linear Embedding
    1. 3.1 Manifolds and Complex Structure
    2. 3.2 Formalizing LLE
    3. 3.3 Applying LLE to the ANES Data
    4. 3.4 A Risk-Averse Workflow
    5. 3.5 Suggestions for Further Reading
  4. 4 Nonlinear Dimension Reduction for Visualization
    1. 4.1 From Linear to Nonlinear Dimension Reduction
    2. 4.2 t-SNE
    3. 4.3 UMAP
    4. 4.4 A Tidy Tangent on UMAP via uwot::umap() and embed
    5. 4.5 Suggestions for Further Reading
  5. 5 Neural Network-Based Approaches
    1. 5.1 A Basic Neural Network Architecture
    2. 5.2 Self-Organizing Maps
      1. 5.2.1 Defining the Steps of the Algorithm
      2. 5.2.2 Applying SOM to the ANES Data
    3. 5.3 Autoencoders
      1. 5.3.1 Applying Autoencoders to the ANES Data
    4. 5.4 Suggestions for Further Reading
  6. 6 Final Thoughts on Dimension Reduction

###### List of Figures

  1. 2.1 Correlation between Features and Trump
  2. 2.2 Correlation between Features and Japan
  3. 2.3 Network Representation of Feature Correlations
  4. 2.2 What is PCA Doing?
  5. 2.6 Linear Regression Approach to Summary
  6. 2.7 Projecting Points onto Components
  7. 2.8 Full PCA Solution for Simulated Points
  8. 2.9 PVE and CPVE for the PCA Solution
  9. 2.10 Scree Plots of Variance Explained (Percentages & Eigenvalues)
  10. 2.11 Biplot of PCA Fit
  11. 2.12 Feature Loadings from PCA Fit
  12. 2.13 PCA Scores by Party Affiliation
  13. 3.1 Manifolds and Complex Structure
  14. 3.3 3D Distributions Across 12 Features
  15. 3.4 Optimal $k$ from Minimal Value for $1-\rho^{2}$
  16. 3.5 LLE Projection Results
  17. 3.6 Comparing LLE to PCA
  18. 3.7 2D Contour Plots Across Eight Original Features
  19. 4.1 Manually Varying Perplexity for t-SNE
  20. 4.2 t-SNE Grid Search Across $\theta$ and Perplexity
  21. 4.3 Comparing UMAP Fits Across Epochs: 500 vs. 20, with $k=5$
  22. 4.4 UMAP Grid Search Across $k$ and Epochs
  23. 4.5 A Tidy UMAP Solution
  24. 5.1 Training Progress for the SOM Fit
  25. 5.2 Searching for 2 Clusters via K-Means
  26. 5.3 Searching for 2 Clusters via Fuzzy C-Means
  27. 5.2.2 Applying SOM to the ANES Data
  28. 5.2.2 Applying SOM to the ANES Data
  29. 5.8 Deep Feature Plots
  30. 5.3.1 Applying Autoencoders to the ANES Data
  31. 5.11 Deep Features 12 and 13 by Party Affiliation

###### List of Tables

  1. 1.1 Question Wording for Feeling Thermometers
  2. 1.2 Feeling Thermometer Summary Statistics: Raw Features
  3. 1.3 Feeling Thermometer Summary Statistics: Imputed Features
  4. 1.4 Comparing Summaries for the Three “Most Missing” Features
  5. 5.1 Confusion Matrix For Deep Neural Network Predicting Party Affiliation

## Acknowledgements

I am indebted to the my colleagues and students in the computational social
science program at the University of Chicago. The UChicago environment has
both stimulated and challenged me at a number of turns, which has undoubtedly
impacted this manuscript for good, whether in conversations with colleagues,
writing code for class, picking up new texts that inspire me, or simply
sitting and thinking in my office.

Further, the series editors Mike Alvarez and Neal Beck have been immensely
encouraging and supportive throughout the entire process. As with my previous
Element, a special thanks is owed to Neal Beck, who has been a constant source
of challenge and support. His friendship and insight have not only made this a
notably stronger manuscript, but have also provided me immense personal
benefit.

I am forever thankful for my precious three daughters and to my wife, Becky,
who continues to deepen her love and support for me as my best friend and
closest companion.

## Replication & Code

All code (DOI: 10.5281/zenodo.4594352) replicating all numerical and visual
results are available at **Link to Cambridge Repository xxxxx.zip file**, and
can be openly accessed online via Github,
https://github.com/pdwaggoner/dimension-reduction-CUP.

##  1 Introduction

The modern era of research is more concerned with collecting and learning from
data than ever before. Whether for influencing decision-making, deepening an
understanding of “process” broadly defined, or some other task, excitement and
preoccupation with data are based in the rapid and massive production of data.
Correspondingly, we are witnessing rapid development of new methods for
efficiently learning from these new data.

Accomplishing this central task of learning from data typically takes shape in
a supervised way, as forecasting and predictions are of high importance. Yet,
supervised machine learning, which is concerned with predicting some labeled
value based on learned patterns from a tuned model, in such a rapidly changing
landscape can be tricky. For example, some labels may not exist (e.g.,
measurement error), the target outcome may be unclear, a priori expectations
on outcome patterns may be nonexistent, or the input features may not be
thoroughly or properly motivated. In such cases, the nature of “supervision”
becomes unclear.

Alternatively, unsupervised machine learning offers a different approach to
the task of learning from data. Unsupervised learning allows data to
essentially speak for itself, where no ground truth or expected outcomes
condition the modeling process, nor does predicting some labeled value.
Rather, unsupervised learning is primarily concerned with uncovering latent,
non-random structure in data. By uncovering this structure, a deeper
understanding of the data, and potentially how it was produced, are possible.
Guided by our main task of learning from data, an unsupervised framework eases
the theoretical burden of researchers relying on theoretical innovation or
some assumed underlying data generating process (which could be incorrect), to
instead defer to the data. In a word, at present we are interested in
exploring data.

Importantly, while data are allowed to speak more freely in an unsupervised
setting than in a constrained supervised setting, unsupervised methods are not
without assumptions as well as choices to be made. Assumptions of different
approaches to reducing dimensionality can result in different views of the
data, as demonstrated and discussed across all techniques covered in this
Element. Similarly, the choices researchers make throughout the process can
also impact the patterns and results that emerge from different versions of
different algorithms. As such, great care should be taken along with
justification of choices made along the way to place results into a proper
exploratory framework. To this end, a theme throughout the Element is
encouraging researchers to try out and compare across several versions of
algorithms, as well as set up tuning grids to search across combinations of
hyperparameters where it makes sense. And at a higher level, it is useful to
remember that unsupervised exploration of data as presented and discussed in
this Element is rooted in an effort to parse signal from noise, which is
extremely common in high dimensional data.

There are two main approaches to unsupervised machine learning: clustering and
dimension reduction. In my previous Element, I covered this first realm of
clustering (Waggoner, 2020). Clustering searches some data space for natural
groupings and patterns, and then seeks to partition the data space in a way
that reflects the underlying structural similarity. Dimension reduction, on
the other hand, though still concerned with recovering structure in data, is
instead interested in doing so by creating a simpler version of the more
complex original version of the full data space. As a result, dimension
reduction represents the structure of the original data in a clearer, more
digestible, and usually simpler way. By combining these two Elements on
clustering and dimension reduction, researchers and practitioners are afforded
a firm understanding of a modern approach to unsupervised machine learning for
addressing a host of social science problems and questions.

Though straightforward to define at a high level, unsupervised machine
learning includes many techniques and algorithms aimed at learning from data
as there are many ways to conceptualize structure. Though some methods are
relatively simple to build and understand, many unsupervised methods can
quickly become complex, both in algorithm construction as well as in
implementation and interpretation.

The goal of this Element, then, is to offer a framework for understanding and
applying dimension reduction in a modern context. In service of this goal, I
detail many algorithms, all of which differ in how they treat and process
data, and thus how they conceptualize structure. Motivated by a common goal
(learn from data) and situated in a common framework (unsupervised machine
learning), the diverse suite of methods covered in this Element offers a
representative picture of the rich diversity of the unsupervised dimension
reduction landscape.

While the task is clear enough, it is important to remember that data are
rarely, if ever, simple and tidy. Though an unsupervised approach to our
central task of learning from data, if properly executed, allows for natural
structure to emerge, the complexities of data make this intuitively-simple
task often complex in application. That is, as data complexity deepens, so too
does the process of method selection, implementation, and interpretation of
the patterns and structure we uncover from the data. The reality of increased
data complexity complicating the modeling process is especially true in a
social science setting, where subjects are often people or institutions
(occupied by people), who are inherently complex and messy. Thus, building
efficient unsupervised learners and then meaningfully interpreting patterns in
a social science context are particularly challenging tasks in modern
applications. Yet, a deeper, unified understanding of dimension reduction can
lead to well-motivated, intentional, and justified decisions made throughout
the modeling process. As a result, social scientists can overcome the hurdles
of data complexity in pursuit of the central goal of learning from data in an
unsupervised way.

###  1.1 Defining the Title

To put meat on the bones of the purpose of this Element and thus the value of
dimension reduction in data analysis, I spend a few paragraphs unpacking the
Element’s title, Modern Dimension Reduction.

First, what are dimensions? Dimensions are variables or features of data.
Mathematically, these are column vectors in some data matrix that, when
increased, also increase the complexities among a set of features existing in
a common data space.111Importantly, data complexity is often defined by both
size and dimensionality. Where size refers to the volume of data in a single
space and often dubbed “big data,” this Element focuses on the dimensionality
part of complex data in light of the scope of methods covered. Yet it is
important to highlight many of these techniques, especially UMAP and
autoencoders, easily adapt to big data settings. So, a high dimensional data
set consists of many features with measured values across some set of
observations (row vectors). This, or any data space can be summarized as
$\mathbf{X}$, which signifies an $N\times P$ data matrix with
$n\in\\{1,\dots,N\\}$ observations (rows) and $p\in\\{1,\dots,P\\}$ features
(columns). Of note, different fields use different terms for $p$, e.g.,
“variables”, “predictors”/“regressors” (in a supervised setting), “inputs,”
“features,” and so on. Throughout this Element, I usually refer to $p$ as
features, because this is a more descriptive name for these components of the
data space, as the measured values record specific features of the
observations, $n$, that exist in $\mathbf{X}$.

Of note, many dimension reduction algorithms can take mixed data types (e.g.,
continuous and categorical). But these algorithms tend to perform best with
continuous features rather than categorical features that have discrete levels
or classes, as categorical features lack variance or nuance that is useful to
understand latent structure. That is, if a feature has values of 0 or 1, then
these are orthogonal, discrete differences. The lack of nuance in the scale
limits the value of dimension reduction as we will see. There may be room to
alter the feature construction (e.g., one-hot encoding). Yet, the value of
dimension reduction as covered in this Element is most beneficial when using
continuous, numeric features. We will come back to this point and also the
role of scaling input features in the following sections.

Now, what is dimension reduction? Dimension reduction is primarily concerned
with taking a high dimensional data space, which typically means $p>4$, and
making a simpler version of it, which typically means $p=2$. And a visual
version of $p=2$ might be a scatterplot, with $X$ and $Y$ axes capturing the
reduced first and second dimensions. We address visualization of dimension
reduction results throughout the Element. But the main idea with dimension
reduction is to embed or represent the high dimensional full data space on a
lower dimensional subspace. Why “represent data” in the first place? Taken as
a whole, the high dimensional, original data space is uninterpretable by
humans such that we cannot readily understand a visual rendering of data in
more than four dimensions.

The process of moving from a high dimensional space to a low dimensional
subspace substantively results in making data being more interpretable,
understandable, and digestible. This value of dimension reduction is why it is
sometimes referred to as “low dimensional embedding,” “mapping,” “lower
dimensional projection,” or “data representation.” The task of dimension,
then, is to learn prominent patterns in the higher dimensional setting, and
then project these patterns onto a lower dimensional subspace. The lower
dimensional space acts as a summary of the full space based on the patterns
naturally existing across all features, $\mathbf{X}$. Importantly, there are
many aspects of this definition of dimension reduction that likely make social
scientists uncomfortable. Namely, by summarizing anything, whether data, a
film, or a baseball game, some information is necessarily discarded or left
out. So too in dimension reduction, which often and intentionally throws away
some data for the benefit of a simpler look at the full, complex original data
space. Using the learned information to move from the higher to the lower
dimensional setting requires a choice of discarding some of the information
deemed non-substantive, which is up to each researcher. For example, principal
components analysis (PCA) is one of the oldest dimension reduction techniques,
and it searches for a lower dimensional version of the data space that
maximizes the total variance. The initial principal component is calculated
for the direction along which the data vary most. Subsequent components are
orthogonal to the preceding components, such that unique variance is captured
in subsequent component calculations. The ultimate choice in most PCA
applications, then, is to decide on the number of components that do a good
(enough) job of characterizing the data. The number of selected components
should be less than $p$, as technically up to $p$ components can be found in
any data space, $\mathbf{X}$. Though unpacked later in the Element, the point
at present is that the choice of selecting a subset of components from among
all calculated components requires the researcher to select some of the data
and discard other parts of the data. This choice is equivalent to saying, “I
am OK moving forward explaining most of the data, but not all of it, for the
benefit of simpler, yet still informative data.”

Finally, though the task of projecting a high dimensional data space onto a
lower dimensional subspace is a common one, there are many ways to
conceptualize patterns in data, and thus many methods for learning these
patterns. For example, we might be interested in reducing dimensionality based
on maximizing similarity across features. Though similarity can be
conceptualized and measured in many ways, e.g., correlation, covariance, or
spatial distance, such a goal would put us in the world of PCA, as previously
discussed. Alternatively, if we suspect an aggregation of many small
neighborhoods of data produce a simplified version of the higher dimensional
space, then we might be in the world of uniform manifold approximation and
projection (UMAP). As discussed later, the goal of UMAP is to first learn the
shape and contours of the manifold (a geometric shape defining a set of data)
underlying the high dimensional data space. Once learned, UMAP translates the
learned manifold to a lower dimensional version of it, based on the learned
distances between observations distributed across the manifold. UMAP has the
added benefit of giving a reproducible solution that efficiently balances
global and local structure defined by the learned manifold. The
reproducibility aspect of UMAP is an immensely powerful extension of another
dimension reduction technique, t-distributed stochastic neighbor embedding
(t-SNE), which is also covered later in the Element. As such, our
conceptualization of similarity, and thus how we treat data during modeling,
will place us in vastly different realms for reducing dimensionality. Despite
the many flavors and differences across techniques, though, the task remains
constant: to learn the structure underlying the high dimensional data, and
then produce a lower dimensional, simpler version of the full, complex data.

As (big) data becomes increasingly complex with demand for sophisticated
computational skills also growing, dimension reduction is a critical skill all
quantitative researchers should know and practice. As one of the core
approaches to unsupervised machine learning, dimension reduction is extremely
helpful for making these widely occurring, complex data spaces more
interpretable and manageable. Whether using dimension reduction as a part of a
broader research program through feature extraction, or on its own to learn
natural patterns in data allowing latent structure to emerge in an intuitive
light, its value is no less diminished.

###  1.2 Running Example: 2019 American National Election Pilot Study

As a running example throughout the Element, I use the 2019 American National
Election Pilot Study data (ANES, 2019).222The American National Election
Studies (www.electionstudies.org). These materials are based on work supported
by the National Science Foundation under grant numbers SES 1444721, 2014-2017,
the University of Michigan, and Stanford University. The ANES is a large, opt-
in survey including complete responses from over 3,000 respondents.

Though many rich features are included in the data (e.g., political
preferences, demographics, etc.), I focus on the battery of feeling
thermometers, which measure respondents’ preferences on a host of topics.
Respondents are asked about how they feel toward some person or concept, and
then asked to record that feeling on a scale from 0 to 100, with 0 being
extremely cold toward the concept and 100 being extremely warm toward the
concept. Though the question wording has evolved over several iterations of
the ANES, the 2019 pilot study used in this Element includes consistent
question wording for all feeling thermometers, How would you rate [topic]?

There are 35 feeling thermometers in the 2019 ANES Pilot survey ranging from
political candidates (e.g., Sanders, Trump, Biden, Harris, etc.) to social
issues and people groups (e.g., transgender, Asians, Muslims, journalists,
etc.) and even institutions (e.g., immigration and customs, NATO, the UN,
etc.) and countries (e.g., Mexico, France, Israel, etc.). The value of these
thermometers for present purposes is many are likely collinear with each other
pointing to common variation across features, and equally contribute to the
complexity of the ANES data space. The assumption here, then, is these feeling
thermometers should project onto some lower dimensional, subspace on the basis
of similarity (however defined). The simpler subspace, then, can be used to
understand natural patterns in the American electorate. Of note, I thoroughly
explore correlation across features at the outset of Section 3, prior to the
treatment of PCA.

In an effort to deepen an understanding of the recovered patterns, I color
data points (respondents) in most visualizations based on stated party
affiliations. Taken with the solutions comprising only a battery of apolitical
feeling thermometers and no feature for party affiliation, my assumption,
which also serves as a naive expectation throughout, is that the structure
underlying responses to these thermometers should take shape in a partisan
way. On average, I expect Democrats to be grouped together and distinct from
Non-Democrats who should also be grouped together. It is worth reiterating
that no party affiliation feature is included in the fit of any algorithms.
Party affiliation is only used to contextualize findings and add clarity to
the recovered latent structure.

Social scientists do not typically work with extremely high dimensional data
in the traditional sense with hundreds or even thousands of dimensions. Yet,
recall that substantively interpreting visual patterns in any dimension
greater than four (or really greater than three in most cases) is a nearly
impossible task. Though not a traditionally “big data” application, the
inclusion of 35-dimensional data in this Element still allows for
demonstration of dimension reduction’s value in social science applications.

For reference, question wording along with the ANES labels for each feeling
thermometer are listed in Table 1.1. The shorthand labels in the second column
will be used throughout the Element.

Table 1.1: Question Wording for Feeling Thermometers Question | Shorthand Label
---|---
How would you rate Donald Trump? | Trump
How would you rate Barack Obama? | Obama
How would you rate Joe Biden? | Biden
How would you rate Elizabeth Warren? | Warren
How would you rate Bernie Sanders? | Sanders
How would you rate Pete Buttigieg? | Buttigieg
How would you rate Kamala Harris? | Harris
How would you rate blacks? | Black
How would you rate whites? | White
How would you rate Hispanics? | Hispanic
How would you rate Asians? | Asian
How would you rate Muslims? | Muslim
How would you rate illegal immigrants? | Illegal
How would you rate immigrants? | Immigrants
How would you rate legal immigrants? | Legal Immigrants
How would you rate journalists? | Journalists
How would you rate NATO? | NATO
How would you rate United Nations (UN)? | UN
How would you rate Immigration and Customs Enforcement (ICE)? | ICE
How would you rate National Rifle Association (NRA)? | NRA
How would you rate China? | China
How would you rate North Korea? | North Korea
How would you rate Mexico? | Mexico
How would you rate Saudi Arabia? | Saudi Arabia
How would you rate Ukraine? | Ukraine
How would you rate Iran? | Iran
How would you rate Great Britain? | Britain
How would you rate Germany? | Germany
How would you rate Japan? | Japan
How would you rate Israel? | Israel
How would you rate France? | France
How would you rate Canada? | Canada
How would you rate Turkey? | Turkey
How would you rate Russia? | Russia
How would you rate Palestine? | Palestine

###  1.3 The Methods

This Element is focused on introducing readers to the intuition of dimension
reduction and its value in applied data analysis and computational modeling
contexts. To this end, along with introduction of hundreds of lines of R code
to guide application and engagement with practical dimension reduction, I
cover six methods.

Principal components analysis (PCA) is a fitting method with which to begin
any treatment of dimension reduction, whether formal or applied. PCA is one of
the earliest approaches to explicitly reducing dimensionality and complexity
of a data space, as opposed to informally reducing dimensionality through two-
dimensional visualization, which has been around for hundreds of years in
various forms. PCA, as tacitly mentioned above, is concerned with finding a
new, lower dimensional subspace based on maximizing the natural variance that
exists across the full set of inputs, $\mathbf{X}$. PCA computes a new set of
components based on a linear combination of weighted features values, which
are often considered as “prototypes” that capture unique variance in higher
dimensions. Observations load onto each component, which are orthogonal to all
other components. The loading of observations onto components results in a new
set of measured values for observations across a new set of features. These
so-called component scores can be extracted and used as new features in
subsequent analyses or can be used for analytical purposes in their own rite.

Next, building on the linear construction of PCA, locally linear embedding
(LLE) is a newer technique based on a similar intuition. LLE linearly combines
weighted feature values to learn a latent structure as in PCA, but it does so
in a neighbor-based way that is more overtly rooted in manifold learning. LLE
searches small neighborhoods in the full data space assuming data are
distributed across a latent manifold, and then calculates weights for
observations existing around the candidate observation, $\forall i\in N$.
Higher weights are given to closer observations, and lower weights are given
to distant observations. Then, at this point, the construction is similar to
PCA, where distance-based weights are multiplied by raw feature values and are
then summed to give a lower dimensional picture of the underlying structure.
LLE learns structure locally and uses the local structure to reproduce the
global structure, whereas PCA is interested in producing a summary of global
structure based on maximal variance.

After PCA and LLE, we abandon linear combinations of weighted feature values,
but build on the manifold aspect introduced with LLE. The focus next is on
explicitly nonlinear approaches to searching some data space for an efficient,
lower dimensional version of the full space. Specifically, we cover
t-distributed stochastic neighbor embedding (t-SNE) and uniform manifold
approximation and projection (UMAP). These approaches are neighbor-based like
LLE, but they differ in how they approximate the lower dimensional version of
the data. Namely, they are both interested in minimizing some cost function
(Kullback-Leibler divergence in t-SNE and cross-entropy in UMAP; discussed
much more later) to reduce dimensionality but based on balancing both local
and global structure learned from the high dimensional space. Then, they
attempt to replicate the learned same structure in a lower dimensional
setting. Though feature extraction is possible with these techniques, they are
much more widely used for visualization in typically two dimensions. I cover
these techniques focusing also on their visual value for learning and
representing the dimension reduction solution.

In the final substantive section, we once again shift gears, and focus on
neural network-based approaches to dimension reduction. These covered
techniques are unsupervised in that we have no expected outcome nor are we
working with labeled data. Instead, we are interested in emulating the way the
human brain learns by firing neurons to give some output based on raw inputs.
This process, as well as it’s connection to computational modeling are
discussed at length in the respective section. The methods are, first, self-
organizing maps (SOM), and then autoencoders. SOM have been around for several
years, and autoencoders are, by comparison, more recent innovations on older
techniques, specifically the restricted Boltzmann machine. SOM are feedforward
neural networks with no hidden layer, and autoencoders are feedforward neural
networks, but with an output layer the size of the input layer. These terms
will be fully unpacked later, but for now, the core idea with autoencoders is
information loss is forced through a process called encoding. The encoded
version of the high dimensional data acts as the “hidden layer.” Then, the
task of the autoencoder is to decode this layer in an attempt to reproduce the
original high dimensional input space, but only on the basis of the simplified
encoded layer. The difference between the output and input layers is called
reconstruction error, and this captures the quality of the autoencoder. It is
important to note the autoencoders can be either shallow or deep, where the
latter version puts us into an explicitly deep learning realm. Thus, we will
cover what it means to move from shallow to deep, and associated benefits and
drawbacks to both.

Also, where it makes sense, important associated innovations or related
concepts will be addressed in an effort to situate readers more firmly in this
world of modern, unsupervised dimension reduction. For example, when covering
LLE, t-SNE, and UMAP, I will place them in the broader subfield of manifold
learning, which is a branch of machine learning dedicated to deriving some
latent manifold, assuming one exists. Or, for example, it is impossible to
understand the final section on SOM and autoencoders without introducing a
basic neural network architecture, which we will do for the task of dimension
reduction.

Of note, there are many other approaches to dimension reduction that are not
covered in this Element, such as multidimensional scaling, ISOMAP, and factor
analysis. Though the main reason for exclusion of other techniques is a
limited amount of space, some of these methods can also suffer from pernicious
problems such as computational inefficiency (ISOMAP compared to UMAP) or a
poor rendering of the high dimensional space (factor analysis compared to deep
autoencoders). In short, though biased by my preferences at some level, I have
made every effort to justify selection of these methods on the basis of
computational efficiency, creativity in dealing with high dimensional data
complexity, generalizability of the solution, as well as linkages across
methods, all of which can be linked back to principal components analysis
(PCA), as we will see at a later point.

###  1.4 The Methods Not Covered

Though we cover a lot of ground in this Element, it is impossible to cover all
approaches to dimension reduction. Of note, I do not address multidimensional
scaling (MDS), other approaches to space scaling (e.g., optimal
classification, NOMINATE, etc.), factor analysis, item response theory (IRT)
models, non-negative matrix factorization (NMF), or Bayesian flavors of these
and many other methods. The reason for this is largely due to insufficient
space, but it is also due to substantive goals. Paired with my recent Element
on clustering (Waggoner, 2020), I am interested in widening the social
scientist’s unsupervised machine learning toolbox. And further, many of the
methods not covered in this Element are either widely used and taught in
quantitative social science applications or are covered eleswhere in other
excellent texts, such as Armstrong et al. (2014).

Also, related to presenting an unsupervised machine learning approach to
dimension reduction, many of the methods not covered in this Element’s
approach dimension reduction from a very different perspective, which is often
more assumption-laden. For example, factor analysis starts with assuming
latent factors are Gaussian to allow for asserting conditional independence
across the factor structure. The goal then is to recover this latent factor
structure that results in the observed input space. That is, factor analysis
assumes there is unobserved structure that causes the observed structure
comprised of input features, e.g.,

$X=b_{1}\text{Component}_{1}.$ (1.1)

In this way, factor analysis is arguably more interested in measurement rather
than explicit dimension reduction, though the line between measurement and
dimension reduction can quickly become gray.

Yet, principal components analysis (PCA) on the other hand stops short of
making any assumptions, distributional or otherwise, related to the “latent
factor structure.” Instead, PCA is interested in learning a reduced version of
the data space comprised of components that are combinations of weighted
features, e.g.,

$\text{Component}=b_{1}X_{1}.$ (1.2)

The goal of PCA, as described in greater depth later in the Element, is to
find the optimal weights (e.g., $b_{1}$ in equation 1.2), which are defined by
maximizing variance on the basis of what we get to observe, which again is the
set of unlabeled input features. Thus, while factor analysis and PCA are both
interested in learning latent structure at some level, these methods consider
and search for this structure in fundamentally different ways, with the former
being a more assumption-heavy approach. Though omission of certain methods in
this Element is largely a practical concern, it remains a second-order goal to
push the reader to think about dimension reduction from a perspective that may
be both new and (hopefully) challenging.

###  1.5 Tidy Programming in R

In this Element, wherever possible, I will leverage the tidy approach to
programming in R (Wickham et al., 2019), which is gaining ground in the social
sciences (Kennedy and Waggoner, 2021). Essentially, the reason for this and at
the core of the tidyverse, is the notion of consistent syntax and stacked
functions. Streamlining key functions across multiple packages with a common
syntax gives way to writing efficient and reproducible code. Of note, the
“stacking” of functions previous mentioned occurs via a key operator: the
pipe, %$>$%, which is read “then…”. The pipe allows for multiple functions to
be passed to each other, which allows for writing large chunks of code to be
run simultaneously, as demonstrated in the following section.

Though readers need not be expert programmers, some familiarity with R, or
other open-source object-oriented languages (e.g., Python or Julia) will
significantly ease reading and implementation of techniques covered
throughout. Where possible, I do my best to annotate and explain code
construction and choices. But I leave it to readers to take what they will
from the code written for this Element.

Readers are also encouraged to load the 2019 ANES Pilot Study data with the
relevant code in the Github repository, which as written requires installation
of the tidyverse and here packages. For the latter package to work, readers
should ensure that the script for running the code in this Element is in a
subfolder called, Data. Otherwise, the function will not work. I highly
recommend inspecting the help documentation for the here package to learn more
about it’s immense value in large scale projects of this sort.

###  1.6 Cleaning the Data

Rather than jump into a technical, though applied work of this sort with clean
data, part of my goal with this Element is to demonstrate an effective
workflow for analysis from start to finish. I strongly encourage readers to
think very carefully about collection and cleaning of data, as much as they
would about modeling and inferences on the back end. Indeed, data, especially
in the social sciences, are rarely clean and complete. Even if complete, data
are rarely in the shape required for specific research projects. In a word,
getting data in good shape is called preprocessing. The first major step in
preprocessing is to clean and organize the data in line with specific project
goals. The second major step, covered in the following subsection, is
addressing missing data. In some circles the second step is called feature
engineering. The code used to preprocess the ANES data used throughout, along
with additional code for techniques not covered in this Element, are included
in the complementary repositories on Github.

There are two main approaches to cleaning data in R: base R operations or the
tidy approach. The base R approach is typically considered manual, e.g.,
selecting features using indexing (e.g., using $\texttt{data\\_set}[\hskip
5.0pt,-3]$ to omit the third feature from the object data_set). A tidy
approach to data management on the other hand, uses consistent syntax with
“human-readable” code to result in a clean, single chunk of code comprised of
many piped functions, which streamlines the research process.

For present purposes, we need to first restrict the space to contain the
features we care about. There are 900 features in the full ANES data set. So,
to make this more manageable and avoid the arduous task of manually selecting
features for each model, the first step is to select() the relevant features
for modeling, which includes the 35 feeling thermometers and an indicator for
party affiliation. This subset of the full feature space is stored in a new
dataset called anes_raw. This is a convenience that will ease plotting both
raw features as well as model objects later in the Element.

With the subset in hand, we need to recode missing values to be a value of NA.
This is a critical step to clarify missing cases in the data, which will be
built upon in the following subsection on feature engineering.

Finally, we take a glimpse() of our cleaned data subset to ensure everything
looks as it should. The glimpse() function from the tidyverse displays feature
names and the first few values, along with displaying the dimensions of the
data matrix, $\mathbf{X}$. We have 3,165 observations along with the 35
feeling thermometers and a feature for party affiliation, which will be used
for coloring the plots throughout the Element. Summary statistics for the
feeling thermometers are presented in Table 1.2.

Table 1.2: Feeling Thermometer Summary Statistics: Raw Features Feature | Complete | Mean | SD | Min | Q2 | Median | Q3 | Max
---|---|---|---|---|---|---|---|---
Trump | 98.9% | 43.87 | 41.43 | 0 | 2 | 38 | 91 | 100
Obama | 99.4% | 53.53 | 37.54 | 0 | 11 | 59 | 91 | 100
Biden | 98.4% | 42.15 | 33.44 | 0 | 7 | 42 | 70 | 100
Warren | 95.8% | 40.51 | 34.36 | 0 | 4 | 41 | 71 | 100
Sanders | 98.6% | 42.15 | 34.88 | 0 | 5 | 42 | 73 | 100
Buttigieg | 90.1% | 38.66 | 30.40 | 0 | 7 | 41 | 61 | 100
Harris | 93.3% | 35.30 | 30.18 | 0 | 4 | 35 | 58 | 100
Black | 99.8% | 70.86 | 23.78 | 0 | 51 | 73 | 91 | 100
White | 99.8% | 70.91 | 22.82 | 0 | 51 | 73 | 90 | 100
Hispanic | 99.9% | 69.00 | 24.20 | 0 | 50 | 71 | 90 | 100
Asian | 99.8% | 69.30 | 23.54 | 0 | 51 | 71 | 90 | 100
Muslim | 99.9% | 51.90 | 29.67 | 0 | 31 | 51 | 75 | 100
Illegal | 99.9% | 42.13 | 32.01 | 0 | 9 | 45 | 68 | 100
Immigrants | 50.1% | 68.33 | 26.71 | 0 | 51 | 72 | 91 | 100
Legal Immigrants | 49.8% | 73.07 | 24.32 | 0 | 55 | 79 | 93 | 100
Journalists | 100% | 49.87 | 31.57 | 0 | 21 | 50 | 77 | 100
NATO | 95.3% | 54.95 | 27.40 | 0 | 41 | 53 | 75 | 100
UN | 97.9% | 50.59 | 30.59 | 0 | 25 | 51 | 75 | 100
ICE | 98.2% | 52.22 | 33.32 | 0 | 23 | 51 | 84 | 100
NRA | 98.1% | 45.87 | 36.54 | 0 | 6 | 49 | 82 | 100
China | 100% | 37.93 | 24.92 | 0 | 17 | 40 | 52 | 100
North Korea | 99.9% | 20.78 | 23.16 | 0 | 2 | 11 | 35 | 100
Mexico | 99.9% | 53.52 | 26.50 | 0 | 38 | 52 | 72 | 100
Saudi Arabia | 99.9% | 32.08 | 24.18 | 0 | 10 | 31 | 50 | 100
Ukraine | 99.8% | 49.08 | 24.37 | 0 | 35 | 50 | 64 | 100
Iran | 99.9% | 26.56 | 24.79 | 0 | 4 | 20 | 47 | 100
Britain | 100% | 69.14 | 23.69 | 0 | 52 | 72 | 89 | 100
Germany | 100% | 60.45 | 24.822 | 0 | 49 | 60 | 80 | 100
Japan | 100% | 64.79 | 23.87 | 0 | 50 | 67 | 84 | 100
Israel | 99.9% | 59.64 | 29.1 | 0 | 43 | 58 | 86 | 100
France | 99.9% | 58.58 | 25.07 | 0 | 47 | 59 | 78 | 100
Canada | 100% | 72.42 | 23.95 | 0 | 56 | 78 | 92 | 100
Turkey | 99.8% | 37.69 | 23.81 | 0 | 18 | 41 | 51 | 100
Russia | 99.9% | 31.43 | 24.98 | 0 | 8 | 30 | 50 | 100
Palestine | 99.7% | 40.13 | 27.44 | 0 | 13 | 47 | 55 | 100

###  1.7 Addressing Missing Data

Most of the algorithms covered in this Element require complete data to give a
reliable lower dimensional version of the data. In fact, some algorithms have
no mechanism to handle missing data, whereas others simply drop these cases if
they are passed to the algorithm. Listwise, or “brute force” deletion can be a
risky strategy and problematic for generalizations drawn from patterns. It is
highly discouraged to simply drop (or delete) data, as data are extremely
valuable, especially in a social science context where it is costly and time
consuming to repeatedly return to a population and draw new samples. As such,
I handle missing data with imputation, which means to replace the instance of
NA with some calculated value that is assumed to be a plausible substitution
for that which the value likely would have been, had it been recorded. Beyond
losing all or most of our data (which alone is sufficient to bypass the
deletion strategy), as previously mentioned a key reason to impute rather than
delete is because many of the algorithms require complete data to run
properly, or at least produce reliable results. To the latter point, some
algorithms may work with incomplete data, opting to essentially ignore these
cases rather than stop. But if the goal is to understand how some complex
higher dimensional data space maps onto a lower dimensional version of the
original complex space, then it does not make sense to reduce dimensionality
of a partially populated data space.

Before imputing data, I recommend first getting to know the data better to
understand the contours of missingness. Such a step will help diagnose the
cause of the missingness, if any, as well as the precise location of
missingness in the full data space. Though omitted from the text, I leverage a
series of visualizations built with the naniar package and present the code on
Github. Inspecting the patterns of missingness allow for more targeted
strategies for dealing with missingness. Once visually explored, but before
getting to the imputation strategy, it is important to consider the cause of
missingness in light of the patterns of missingness. Rubin (1976) defined
three major causes: missing not at random (MNAR; a systematic issue in
recording or measuring data resulted in missingness), missing at random (MAR;
missingness is random, but is not an equal probability for each observation),
and missing completely at random (MCAR; missingness is equally random for all
observations). Though often difficult to defend one of these underlying causes
of missingness in the data, researchers should take care in at least thinking
about this underlying driver before proceeding to imputation or any strategy
aimed at dealing with the missingness.

Assuming MCAR, I proceed with imputation by taking advantage of the recipes
package, which is has many powerful functions for data preprocessing of this
sort. Specifically, I impute missing values using the k-Nearest Neighbors
algorithm, or kNN. kNN is a simple supervised learning algorithm that, as
leveraged in the current case, imputes based on a subset of observations
surrounding a data point. First, kNN defines a small neighborhood of size $k$
around a candidate case with a missing value. Averages are then taken based on
the values of those existing in the smaller neighborhood. The neighborhood
average is the imputed value for the missing case. Put differently, kNN is
used to first reduce the search space. Then, from among the reduced group of
observations, the average of the non-missing cases in the neighborhood serve
as the new, imputed value for the missing case. The procedure follows these
steps for all missing cases across all features. This approach is considered
multiple imputation in that there is a statistical procedure resulting in a
computed value, rather than the more common single imputation approach that
might be on the basis of taking the same value from the most recent complete
case along each feature. kNN for imputation has the added benefit of deciding
on which features to create the neighborhood. For example, should we define a
neighborhood on the basis of people-based feeling thermometers? Or
institution-based thermometers? Or all features? For this application, I base
imputation on all features. The summary statistics for the full, imputed ANES
data set based on our kNN recipe is presented in Table 1.3. We can see from
the table, first and foremost, that all completion rates are 100%, suggesting
we have successfully imputed the missing cases.

Table 1.3: Feeling Thermometer Summary Statistics: Imputed Features Feature | Complete | Mean | SD | Min | Q2 | Median | Q3 | Max
---|---|---|---|---|---|---|---|---
Trump | 100% | 43.72 | 41.31 | 0 | 2.0 | 37.0 | 90 | 100
Obama | 100% | 53.48 | 37.48 | 0 | 11.0 | 59.0 | 91 | 100
Biden | 100% | 42.26 | 33.26 | 0 | 8.0 | 43.0 | 70 | 100
Warren | 100% | 40.74 | 33.95 | 0 | 4.0 | 41.0 | 70 | 100
Sanders | 100% | 42.26 | 34.74 | 0 | 5.0 | 43.0 | 73 | 100
Buttigieg | 100% | 38.97 | 29.67 | 0 | 9.0 | 42.0 | 60 | 100
Harris | 100% | 35.57 | 29.68 | 0 | 5.0 | 36.0 | 57 | 100
Black | 100% | 70.86 | 23.77 | 0 | 51.0 | 73.0 | 91 | 100
White | 100% | 70.91 | 22.80 | 0 | 51.0 | 73.0 | 90 | 100
Hispanic | 100% | 69.00 | 24.18 | 0 | 50.0 | 71.0 | 90 | 100
Asian | 100% | 69.30 | 23.52 | 0 | 51.0 | 71.0 | 90 | 100
Muslim | 100% | 51.90 | 29.66 | 0 | 31.0 | 51.0 | 75 | 100
Illegal | 100% | 42.14 | 32.00 | 0 | 9.0 | 45.0 | 68 | 100
Immigrants | 100% | 67.52 | 23.69 | 0 | 51.0 | 70.0 | 87 | 100
Legal Immigrants | 100% | 72.67 | 21.43 | 0 | 58.4 | 76.2 | 90 | 100
Journalists | 100% | 49.87 | 31.57 | 0 | 21.0 | 50.0 | 77 | 100
NATO | 100% | 54.85 | 27.09 | 0 | 41.0 | 53.0 | 75 | 100
UN | 100% | 50.63 | 30.40 | 0 | 25.6 | 51.0 | 75 | 100
ICE | 100% | 52.06 | 33.12 | 0 | 24.0 | 51.0 | 84 | 100
NRA | 100% | 45.72 | 36.34 | 0 | 6.0 | 49.0 | 81 | 100
China | 100% | 37.93 | 24.92 | 0 | 17.0 | 40.0 | 52 | 100
North Korea | 100% | 20.78 | 23.15 | 0 | 2.0 | 11.0 | 35 | 100
Mexico | 100% | 53.52 | 26.50 | 0 | 38.0 | 52.0 | 72 | 100
Saudi Arabia | 100% | 32.09 | 24.17 | 0 | 10.0 | 31.0 | 50 | 100
Ukraine | 100% | 49.07 | 24.35 | 0 | 35.0 | 50.0 | 64 | 100
Iran | 100% | 26.55 | 24.78 | 0 | 4.0 | 20.0 | 47 | 100
Britain | 100% | 69.14 | 23.69 | 0 | 52.0 | 72.0 | 89 | 100
Germany | 100% | 60.45 | 24.82 | 0 | 49.0 | 60.0 | 80 | 100
Japan | 100% | 64.79 | 23.87 | 0 | 50.0 | 67.0 | 84 | 100
Israel | 100% | 59.63 | 29.09 | 0 | 43.0 | 58.0 | 86 | 100
France | 100% | 58.58 | 25.07 | 0 | 47.0 | 59.0 | 78 | 100
Canada | 100% | 72.42 | 23.95 | 0 | 56.0 | 78.0 | 92 | 100
Turkey | 100% | 37.68 | 23.79 | 0 | 18.0 | 41.0 | 51 | 100
Russia | 100% | 31.42 | 24.98 | 0 | 8.0 | 30.0 | 50 | 100
Palestine | 100% | 40.13 | 27.42 | 0 | 13.0 | 47.0 | 55 | 100
Democrat | 100% | 0.42 | 0.49 | 0 | 0.0 | 0.0 | 1 | 1

To inspect the quality of our solution and to avoid repeatedly comparing back
and forth with the raw ANES summary statistics in Table 1.2, we can take a
closer, individual look at some of the particularly problematic features. For
these three features with the most missing values, we are looking for
substantively similar summary statistics for the raw version of the feature
compared to the imputed version. Take a look at the comparisons for the
Immigrants, Legal Immigrants, and Buttigieg features presented in Table 1.4.

Table 1.4: Comparing Summaries for the Three “Most Missing” Features Feature | Complete | Mean | SD | Min | Q2 | Median | Q3 | Max
---|---|---|---|---|---|---|---|---
Immigrants (raw) | 50.1% | 68.33 | 26.71 | 0 | 51 | 72 | 91 | 100
Immigrants (imputed) | 100% | 67.52 | 23.69 | 0 | 51.0 | 70.0 | 87 | 100
Legal Immigrants (raw) | 49.8% | 73.07 | 24.32 | 0 | 55 | 79 | 93 | 100
Legal Immigrants (imputed) | 100% | 72.67 | 21.43 | 0 | 58.4 | 76.2 | 90 | 100
Buttigieg (raw) | 90.1% | 38.66 | 30.40 | 0 | 7 | 41 | 61 | 100
Buttigieg (imputed) | 100% | 38.97 | 29.67 | 0 | 9.0 | 42.0 | 60 | 100

In general, when the missingness is greater as for the Immigrants and Legal
Immigrants features, we will do a slightly worse job of matching the summaries
for the raw (non-imputed) version as we simply have less data to learn from in
the full data space. Yet, even still, with completion rates of 50.1% and 49.8%
for Immigrants and Legal Immigrants respectively, which equates to about half
of the values missing, the summary statistics are remarkably close to the
original, raw values. This suggests that the other respondents in the
neighborhood of the missing cases were quite helpful (and similar along other
dimensions) in imputing these missing values. Regarding the third feature,
Buttigieg, the summaries were much more closely aligned with the original, raw
summaries. Substantively, the similarity across patterns of summaries suggests
that respondents are revealing relatively consistent signals about preferences
across these feeling thermometers, regardless of the issues. The redundancy of
signaled information across the feeling thermometers allows for an ideal
scenario to leverage dimension reduction to rid our data space of
redundancies, and instead retain the most unique, useful information and
patterns underlying the data space.

##  2 A Classic Approach to Dimension Reduction

The main idea with dimension reduction is to reduce complexity (that is,
dimensionality for present purposes) of a data space to create a lower
dimensional representation of that original space. Such a task makes data more
accessible, patterns more intuitive, and thereby eases to task of detecting
and interpreting natural structure.

Perhaps the most commonly taught dimension reduction technique is principal
components analysis (PCA), which is a linear approach to dimension reduction
that reduces complexity by maximizing variance across the data space. The
starting place with PCA is assuming variance is the best way to conceptualize
structure in a data space. As such, variance and specifically shared variance
will be a key theme in covering this classic approach to dimension reduction.

Applying PCA requires some information to be willfully thrown out. This might
make some feel uncomfortable, but our main goal of dimension reduction in the
first place is to simplify a complex higher dimensional space. This could be
as simple as moving from 2 dimensions to 1, or as complex as moving from 5000
dimensions to 2. Regardless of the complexity of the problem, by reducing
dimensionality of the data space, information is by definition intentionally
being discarded as we decide what the lower dimensional configuration should
look like. As such, when dimension reduction is the goal, this “downside” of
discarding information for the sake of a simpler representation is in reality
usually not considered a downside at all, because we are assuming the
simplified space retains enough of the information from the higher dimensional
space to allow for an understanding of the structure of the space. Put
differently, and more bluntly, dimension reduction in these terms focuses on
the most interesting parts of the data, and gets rid of the less interesting
parts. In so doing, we are able draw conclusions about latent structure. In
terms of PCA, then, this structure is defined by shared variance across all
features.

###  2.1 Why PCA?

Before diving into PCA, why does it make sense to leverage PCA for dimension
reduction? As referenced in different ways to this point, PCA makes a high
dimensional space less complex, and thus more interpretable. With this overt
benefit, come several additional benefits such as visualization, understanding
patterns underlying the data, and “feature extraction” for creating new
features to be used for other tasks, e.g., prediction.

Further, PCA helps with diagnosing and bypassing limitations caused by
multicollinearity across the feature space. As hinted at in the previous
section, the feeling thermometer space in the ANES data is likely
characterized by a sizeable amount of overlap in responses to the battery of
issues. For example, feelings toward Joe Biden, Barack Obama, and Elizabeth
Warren are likely picking up on commonality in preferences of respondents, as
these people are all Democrats who ran for the U.S. presidency, and who have
commanded active media attention. Similarly, feelings toward the U.S.
Immigration and Customs Enforcement Agency (ICE) and feelings toward
immigration likely correlate as well, either negatively or positively, as
these are addressing different aspects of a common underlying concept:
immigration. The expectation with features that correlate then, is that they
also explain similar variance in preferences of people who responded to the
ANES survey. Substantively, many of the features in the full feeling
thermometer space might be more clearly and parsimoniously characterized when
we focus on uncovering the underlying, shared patterns (e.g., immigration),
which is ultimately recorded in a reduced version of the data space that is
picking up on these shared aspects of respondents’ preferences.

Further, and closely related, dimension reduction helps guard against the
curse of dimensionality, which asserts that as the dimensionality of the space
increases, true similarities and differences across the dimensions/features
becomes less clear. To guard against obscuring real similarities and real
differences that naturally exist in the data in higher dimensions, dimension
reduction, which simplifies this complex space by removing dimensions, is an
extremely useful tool.

To make these ideas, and thus the benefits of dimension reduction come alive,
I shift to show correlations across all feeling thermometer features to help
us clarify this assumption of likely shared variance and correlation across
features. Importantly, correlation and variance are distinct statistical
concepts, which are unpacked in the formal definition of PCA below. But to
this point, they are used interchangeably to draw out the substantive value of
PCA for dimension reduction. Another way to put this is to begin by diagnosing
the space to get a clue as to whether structure and similarity naturally exist
in the feeling thermometer space.

In this demonstration, several view of correlations across all feeling
thermometers are provided. First, I compare all features to feelings toward
Donald Trump. The justification for doing so is purely descriptive, where
Trump (being a Republican president) and frequently in the news, is a
political figure allowing us to check some base expectations. For example, we
might expect feelings toward Barack Obama (presidential predecessor and member
of the opposite party) would be strongly negatively correlated with feelings
toward Trump. To do so, we load the cleaned version of the ANES data from the
previous section, and then work primarily with the corrr package from the
tidyverse. Importantly, as corrr is apart of the tidyverse, we can directly
pipe plotting functions giving a clear rendering of the correlations in Figure
2.1.

Examine the correlations between all features and feelings toward Trump in
Figure 2.1. Indeed, as expected, several features correspond with an
intuitive, base set of expectations. For example, feelings toward Barack Obama
are indeed most strongly and negatively correlated with feelings toward Donald
Trump. And feelings toward the National Rifle Association (NRA) most strongly
and positively correlated with feelings toward Trump.

Figure 2.1: Correlation between Features and Trump

Whereas feelings toward Donald Trump took shape in a largely intuitive way
given the salience of Trump’s presidency and the media attention he commands,
we turn now to a different case, and check correlations in another light by
exploring across all features in relation to feelings toward Japan, where
expectations of correlation patterns are perhaps less obvious. See the results
in Figure 2.2.

Figure 2.2: Correlation between Features and Japan

Notably, correlations between all features and Japan are relatively strong and
positive, with the exception of feelings toward Trump. Though interesting and
confusing, possible reasons for exactly why this is the case are beyond the
scope of present purposes. Rather, I am interested in demonstrating that there
are clearly strong correlations across the data space, some intuitive (Trump)
and others less so (Japan).

We continue this exploration of correlations naturally in the feature space,
which will deepen an answer to our motivating question on why it is useful to
pursue PCA to simplify this space. Rather than explore correlations between a
single feature and all others, we can instead view a network representation of
the correlations that naturally exist in feature space, which provides a
fuller, more nuanced picture of degrees and directions of correlations across
all features. To do so, we still rely on core tidyverse packages, corrr and
ggplot2. But for this case, we leverage a different function, network_plot(),
which gives the network. Further, we use the amerika package to color the
network (Waggoner, 2019), as well as many of the visualizations used
throughout the Element.

Figure 2.3: Network Representation of Feature Correlations

A few useful trends emerge from the network configuration of correlations in
the full features space shown in Figure 2.3. First, there is strong
correlation (darker shades) across much of the plot, implying strong
correlations across the full high dimensional data space. Related, a clear
structure seems to characterize this space in a way we might expect, whether
an “expert” or not in American politics. That is, feelings toward countries
are largely grouped together in the upper part of the plot, feelings toward
issues are largely grouped together in the lower right of the plot, and
feelings toward people tend to be grouped together on average in the lower
left of the plot in Figure 2.3. These groupings emerge from the construction
of the network_plot() function, which groups based on the strength of
correlations within a subgroup of features. The goal, really, is to push these
groups apart from one another, again on the basis of natural structure, to
essentially exaggerate differences between features on the basis of the
strength of the correlations. This plot is similar to the concept of
modularity in graph theory, where stronger ties within a module/cluster as
well as sparser connections between modules or clusters implies latent
structure in a common space. Similarly, some latent structure based on clear
collinearity across features is present in this space. The results across
Figures 2.1, 2.2, and 2.3, then, offer sound motivation to move forward with
dimension reduction.

###  2.2 What is PCA Doing?

We have a clear sense that the feeling thermometer feature space is indeed
highly correlated. With this information, we might conclude that it makes
sense to progress with PCA to draw out and isolate dimensions of greatest
variance underlying the data. Yet, before formalizing PCA, we must first
address what it is doing, or precisely how it handles the task of dimension
reduction. I prefer to start in words, and then use equations to clarify.

When approaching a task of summarizing data for the purpose of making it
simpler, which is at the heart of dimension reduction, there are many ways to
go about this. For example, we might overlay some line to summarize a bunch of
observations. There are many possible lines that could exist as a summary,
with some better than others. Consider the hypothetical case of five data
points in Figure 2.5, along with some candidate summary options in Figure 2.5.

Figure 2.4: Five Data Points

Figure 2.5: Summary Options

Such an incremental approach to honing in on an optimal summary of data shown
in Figure 2.5 not only inefficient, but also makes it impossible to hone in on
an optimal option absent a definition of optimality. Simple linear regression
handles this problem in a parametric way, by estimating parameters (an
intercept and a slope) that minimize the sum of squared residuals, on the
basis of which we can represent this summary with a line through the data,
with some uncertainty of course, given the process of parameter estimation.
See this approach in Figure 2.6.

Figure 2.6: Linear Regression Approach to Summary

Rather, PCA, which is technically a special case of linear regression, but
with no estimated parameters and thus no intercept, approaches this problem
from a different perspective with a different goal in mind. In PCA, the focus
is on summarizing a higher dimensional data space by focusing on maximizing
the total variance of that space. As variance is the focus, thereby giving us
a definition of optimality, PCA is initialized to look for the direction in
the data along with the full data vary most. Once that direction is found, PCA
computes and places a summary line called a principal component, that
summarizes the direction of maximal variance. We can call this first principal
component $C_{1}$.

If we stopped with $C_{1}$, we would be explaining a good amount of variance
in the data in most cases, though not the total amount by definition.
Variance, especially in high dimensional contexts, can be complex and proceed
in many directions. Thus, once $C_{1}$ is found, PCA proceeds to search for
the next direction in the data along which the second-most amount of variation
exists. We can call this second principal component $C_{2}$, which now
summarizes the unique, remaining variance after accounting for $C_{1}$. This
process is continued until we have explained the full data space, such that
$C_{*}=p$. Sticking with the number of components the size of the
dimensionality of the full data space is equivalent to saying we have
summarized the entire amount of variance in original setting. At this point,
it would not make any sense to continue with PCA, as we would be back in the
high dimensional setting, where $C_{*}=p$. Instead, our task is to home in on
some subset of the components, $C<p$, to accomplish our goal of simplifying
the high dimensional, complex data space.

As the algorithm searches for and places new prototype components, $C$, to
summarize unique variance in the full data space, uniqueness is defined by a
requirement of each component being orthogonal to all others. This means we
can only define a new principal component if it is explaining completely
unique, previously unexplained variation in the data.

Importantly, as we are operation in some $d$-dimensional space, placement of
components is defined both by feature values and observation placements in the
high dimensional setting. That is, part of the component placement is
constrained by the projection of individual points onto the component. In so
doing, we are left with component scores across each component. Returning to
our simple example in Figure 2.5, we can see what this process looks like for
two components in Figure 2.7.

Figure 2.7: Projecting Points onto Components

The left plot in Figure 2.7 is showing the first principal component (first
dimension), which is summarizing the maximal variance in the data, which in
the case is a pattern that extends diagonally from the upper left of the plot
to the lower right of the plot. Individual points are projected onto the
component, giving these points new measured values in the first dimension of
the new lower dimensional space. Think of these component scores the same as
the measured values for any feature like self-reported political ideology.
Then, the right plot in Figure 2.7 shows the second component, which is
orthogonal to the first. Similarly, points are projected onto the component,
and these scores are the new “measured” values for the second dimension. If a
researcher stopped at this point, then $C1$ and $C2$ could be used as input
features for some predictive modeling task.

Thus, upon finding these components, we would get the full solution shown in
Figure 2.8. Figure 2.8 more clearly shows the uniqueness aspect of the PCA
solution, where each component is orthogonal to the other.

Figure 2.8: Full PCA Solution for Simulated Points

###  2.3 Formalizing PCA

With the previous substantive discussion in mind, we can formalize these ideas
with equations. Recall I mentioned in the previous subsection that PCA is a
special case of linear regression, just without estimated parameters or an
intercept. We can see that in the base construction of PCA for the first
component, $C_{1}$,333I use $\beta$ to emphasize the link to linear
regression, given the wide spread familiarity with regression. The weights are
essentially operating in very similar ways as an estimated $\beta$ coefficient
in linear regression.

$C_{1}=\beta_{11}X_{1}+\beta_{21}X_{2}+\ \dots\ +\beta_{p1}X_{p}$ (2.1)

where, $\beta_{11}$ is a weight for the first component on the raw feature
$p\in\\{1,\dots,P\\}$. From this form, we can see the unique contribution of
each feature to the calculation of each component, $C_{*}$. These weights are
typically called loadings, which captures the idea of each raw feature, $p$,
loading differently onto each calculated component, $C_{*}$. The relation to
the simple linear regression should be clear from Equation 2.1, where the
component, which recall is also a summary line that passes through the data
and the origin, is a function of up to $p$ unique contributions.

Importantly, with PCA and all dimension reduction techniques covered in this
Element, standardization of the features is critical, given unique raw feature
values can take over variance of other features if these features are on
different scales. For example, if there were two features weight in pounds and
income in U.S. dollars, the vastly different scale of the features would
result in an imbalance of contribution to the calculation of the component
solely on the basis of differences across their scales, rather than unique
variance in relation to all other features. To account for this issue,
standardization, which is simply defined by dividing each feature by it’s
standard deviation, guards against this threat. The result is all features are
allowed to be directly compared in a scale-free way.

Then, for each subsequent component, we simply update the index notation
(e.g., $\beta_{11}=$ feature 1, component 1), and find the next component,
$C_{*}$, orthogonal to and thus uncorrelated with the preceding components.
Such a strategy allows us to continue to pick up unique, left-over variance
with our solution,

$C_{*}=\beta_{1*}X_{1}+\beta_{2*}X_{2}+\ \dots\ +\beta_{p*}X_{p}$ (2.2)

James et al. (2013) offer a simple framing of the problem of finding optimal
$\beta$ values for each component, $C_{*}$, to maximize the sample variance
across the observations, $i\in\\{1,\dots,n\\}$, and features,
$p\in\\{1,\dots,P\\}$,

$\frac{1}{n}\sum_{i=1}^{n}\left(\sum_{p=1}^{P}\beta_{p1}X_{ip}\right)^{2}.$
(2.3)

Viewed as an optimization problem, finding optimal $\beta$ values in the
loading vector, $\beta_{*}=\\{\beta_{1*},\dots,\beta_{p*}\\}$, can be solved
using many techniques such as singular value decomposition (SVD). Overly
simplified and in words, SVD is generally comprised of three steps. First,
project the observations on the component, and store the coordinates. Second,
calculate the distance from each point to the origin, which always has
cartesian coordinates, $(0,0)$. Third, square each of these values and add
them together. This series of steps gives the eigenvalue (EV) for each
principal component. Values are squared, because $\sqrt{EV}$ is equal to the
singular value, which is involved in the decomposition of $\mathbf{X}$. For a
more thorough treatment of SVD, see chapter 14.5 in Friedman, Hastie and
Tibshirani (2001).

A final step in a PCA solution is to decide on the number of components to
retain, which is the step referenced a few times to this point deciding “how
much information to throw out and how much to keep.” Though there is no formal
rule for determining this, there are a number of descriptive techniques that
can help. But before getting to these, a final definition that is integral to
understanding PCA is the proportion of variance explained (PVE). The PVE is a
normalized (to equal 1 across all summed components) value that gives an
indication of each component’s contribution to the full PCA solution. Again,
though no clear rule exists for evaluating these as it relates to determining
the optimal number of components to retain, it is reasonable to suggest a
total PVE of around 75-80% is a fair base line, as this could include either a
single component that is doing the bulk of the explanation, or it could
include several components that are contributing to a simpler version of the
high dimensional space.

###  2.4 Applying PCA to the ANES Data

With a clear idea of why PCA is useful, what PCA is doing, and how PCA works
to reduce dimensionality, we conclude this section with an application of PCA
using the 2019 cleaned ANES data. Once applied, we will discuss the output and
the several options for evaluation and interpretation.

There is a remarkably small amount of code needed to fit a PCA model in R. We
will be using the prcomp() function from base R, given the long-standing
status of PCA in applied statistics and statistical computing. Some other
packages have PCA functions such as the FactoMineR or adea4 packages. Yet, in
practice it is much more common to use prcomp() to fit a PCA model given it’s
computational efficiency and simplicity. As such, this is where we start as
well.

We will work with a new package for excellent, simple plotting options, which
is built upon the tidyverse’s ggplot2 package. With the data loaded, we fit
the PCA model on all feeling thermometers in the cleaned dataset, withholding
the party feature, democrat in the 36th column of the data matrix, to use for
visualization of PCA scores.

When we summarize the model output, the “Importance of Components” is
returned. The PVE, which is the second row of values, tells us the proportion
of variance explained by each calculated component. Recall, the PVE tells us
how much of the unique variance is explained by each of the PCs. As expected,
the PVE is decreasing as we move from left to right as the requirement for
defining uniqueness in PCA is defined by subsequent components being
uncorrelated and orthogonal to all preceding components. Thus, we are left
with progressively less variance to be explained as we continue to find and
calculate components. Related, note that the solutions returns 35 components
(denoted by $PC1$, $PC2$, and so on). This is the case, because as previously
mentioned, when we explain all of the variance in a data space, we are back in
the high dimensional setting, which by definition is fully explaining itself
on the basis of the inclusion of all features. The task of PCA, then, is to
hone in on a reduced version of the full space on the basis of explained
variance. The PVE, naturally, helps us out with this task.

Next, and related, we can see the cumulative PVE values (CPVE). These values
are the inverse of the PVE, where they can be progressively summed, and will
eventually equal 1, once summed across all components. For example, we see a
clear jump with no component ($CPVE=0$) to the $PC1$, which has a
$CPVE=0.3599$. Then, when we proceed to $PC2$, we get an increase in PVE of
$0.1657$, as CPVE is at 0.5256 when accounting for PC2, minus the PVE 0.3599
based only on PC1, gives an increase of 0.1657 moving from PC1 to PC2. Summing
PVE for each subsequent PC, we get a CPVE of 1.000 at $PC35$, meaning we have
now explained the full data space, such that $C=P$.

Finally, the model output returns the standard deviation (first row). Recall,
in statistics the standard deviation is a measure of spread and is defined as
the square root of the variance. And recall also, we previously noted that in
PCA the variance is defined by the eigenvalues across the components. And
finally, recall that we said the singular value is simply the square root of
the eigenvalues. Thus, we interpret the standard deviation from PCA output as
the square root of the eigenvalues computed for each principal component. The
decrease in standard deviation values as we move from the first to the final
principal component, thus, make sense, as we are left with progressively less
variance to explain once we reach $C=P$.

We can also call the loadings, which are feature contributions to each
principal component, by calling pca_fit$rotation. The output is omitted due to
its size. Yet, a more effective method for evaluating PCA output is visually.

To do so, we start by visually evaluating the structure of the space, which
builds on the previous numeric exploration of the PCA output. We will manually
calculate the PVE and CPVE, and create two ggplots for each respectively.
These plotted values over each component are included in Figure 2.9. Running
the respective code will first make the calculations and store the values
accordingly. Then, the subsequent code will produce both plots side-by-side,
with labels according to the component number, ranging from 1 to 35 and
created using the ggrepel package.

Figure 2.9: PVE and CPVE for the PCA Solution

Indeed, in Figure 2.9, we get visual corroboration with the numeric output,
that the first few principal component explain the majority of the variance.
Even the first two components explain over half of the data, suggesting a
large amount of correlation in the full feature space.

Another popular approach is to use a scree plot to evaluate the dimensionality
of a space. For this approach, the factoextra package includes options for
scree plots for either percentage of variance explained, or eigenvalues for
each principal component. For either version, which give identical information
as each are capturing differences in variance by component, simply call the
fviz_screeplot() function, and specify the proper choice, either “variance” or
“eigenvalue.” See these results in Figure 2.10.

Figure 2.10: Scree Plots of Variance Explained (Percentages & Eigenvalues)

Figure 2.10 shows the first few dimensions/components seem to be explaining
the bulk of the variance.

Though we have options to explore variance explained by the PCA fit, we need
to use this information to determine how many dimensions to retain in our
simplified version of the data space. Some make suggestions based on total
variance explained as previously discussed, and others suggest components with
eigenvalues greater than 1. Still others suggest looking for the “knee” or
“elbow” of the scree plot to make a get decision. If all of these approaches
sound murky, that is because they are. There is no formal guidance on the
optimal number of components to retain to define the lower dimensional space.
The best we can do is inspect the results in several ways, as we have done to
this point, and then make a decision. Thus, across all of these suggestions, I
would conclude that likely 4 dimensions characterize the space well, as we see
a drop off in PVE, percentage, and eigenvalues around four and five
components. Given that anything greater than 4 dimensions is virtually
uninterpretable by our brains as discussed earlier in the Element, four
dimensions seems like a reasonable cut off. Even still, the first two
dimensions still capture a large amount of variance, which will be useful for
plotting component scores at the conclusion of this section.

Before we get to scores, though, two other useful ways to visually assess a
PCA model, are a biplot and a plot of the feature loadings. We will walk
through both using the factoextra package again for clean, simple code.

First, the biplot of the PCA fit in Figure 2.11 plots the scores in two
dimensions, where dimension 1 is along the x-axis and the PVE is in
parentheses, and dimension 2 is along the y-axis. The points are the component
scores, which recall in two-dimensional space is the coordinates for the
projection of points onto both principal components. The arrows in a biplot
show the connection between features and each dimension. To create the biplot
in Figure 2.11, run the respective code for Section 3.

Figure 2.11: Biplot of PCA Fit

The dashed lines at $(0,0)$ in Figure 2.11 are for reference only. In Figure
2.11, a clear pattern emerges building on the earlier base expectations
regarding Trump and other concepts often associated with Trump (e.g., Trump,
NRA, ICE, White, Russia, Israel) are characterizing the second component. The
first component, though, is more diffusely characterized by most of the other
feeling thermometers, though extremely closely by feelings toward Iran,
Palestine, and Muslims. Thus, from this simple rendering of the PCA solution,
we can start to get a hint of similarities across features, and how these can
be more simply represented in a lower dimensional setting.

Next, consider the other approach to visually interpreting PCA output by
plotting the loadings Figure 2.12.

Figure 2.12: Feature Loadings from PCA Fit

To help us interpret the output from Figure 2.12, it is useful to point out
that a constraint is placed on the search for the optimal loading values in
Equation 2.3. That is, we are interested in maximizing sample variance across
the data space, but subject to a normalization factor,

$\sum_{p=1}^{P}\beta_{p1}^{2}=1.$ (2.4)

That is, the sum of the square loadings must add up to one. Though James et
al. (2013, 376) and others emphasize this normalization constraint in the
context of preventing “arbitrarily large variance,” another useful result from
normalization of this sort is in the ease of interpretation. Restricting the
size of the coefficients in such a way allows us to effectively interpret the
loadings downstream as we might a correlation coefficient, especially because
we also get negative and positive loadings where features load positively or
negatively onto different components.

As such, interpreting Figure 2.12 showing the feature loadings on each of the
first two dimensions, we start by inspecting the tip of the arrow. At the tip
of the arrow, the contribution of the feature to the component’s calculation
is either negative if it is to the left of (below) 0.0, or positive if to the
right of (above) 0.0. Shorter arrows, then, reflect less correlation with, or
contribution to, the dimension(s). Longer arrows reflect greater correlation
with, or contribution to, the dimension(s). For example, in Figure 2.12, the
Israel feature negatively loads onto the second component to a degree of about
$-0.75$. All input features can be interpreted accordingly with a feature
loadings plot. And of note, readers can double check the loading values by
calling them directly from the PCA fit via pca_fit$rotation.

We will conclude the PCA section with a return to our guiding goal from the
outset of the Element, where we are interested in exploring whether the latent
structure in the feeling thermometer space varies along a partisan dimension.
To do so, I generate a custom view of the solution along the first two
dimensions by plotting the scores against each other and then coloring the
points by party affiliation. Recall, the pid7 party affiliation feature in the
original ANES data included eight levels, where $1=\text{Strong Democrat}$,
$2=\text{Moderate Democrat}$, $3=\text{Lean Democrat}$,
$4=\text{Independent}$, $5=\text{Lean Republican}$, $6=\text{Moderate
Republican}$, $7=\text{Strong Republican}$, and $8=\text{NA}$. To simplify the
plots, I recoded this feature to be dichotomous, where $1=\text{Democrat}$,
and $0=\text{Non-Democrat}$. The goal with this step is to avoid discarding
data or information (e.g., dropping all non-Democrats or non-Republicans),
while instead grouping those who identify at any level with the Democratic
party ($\texttt{pid7}==1:3$), or do not ($\texttt{pid7}==4:7$). But
ultimately, the chief benefit here is to clarify the visual patterns from the
algorithmic output, which is a strategy adopted throughout the Element. As
such, the PCA results with color according to party affiliation are shown in
Figure 2.13. Ellipses are loosely drawn around each party group for
descriptive value only.

Figure 2.13: PCA Scores by Party Affiliation

In Figure 2.13 a clear partisan pattern emerges based only on survey responses
to the batter of feeling thermometers. And perhaps more strikingly, recall the
first two dimensions from the PCA solution account for just over half of the
PVE, 52.6%. Thus, our ability to pick up on a clear partisan distinction based
on only half of the variance in the data suggests that, while responses to
these feeling thermometers are not overtly on partisan terms (i.e., “as a
Republican, rate your feelings on X”), there is a pronounced undercurrent of
partisanship in these survey responses; a pattern we will reference back to
throughout the remainder of the Element.

Perhaps the biggest weakness of PCA is the ambiguity surrounding the optimal
number of reduced features to retain, as there is no formal guidance for
making this choice. As a result, the benefits of dimension reduction via PCA
become less clear. Related, in smaller dimensional contexts (e.g., 5 or 6),
the value of PCA drastically diminishes as well. Recall, we can calculate up
to $p$ principal components, such that the total variance can be explained by
the PCA solution. When all of the variance is explained, we are back in the
high dimensional context. Here again, the value of dimension reduction via PCA
is substantially less clear. Thus, PCA is most useful when there is both high
correlation across the full input space such that a reduced set of features
contributes to the goal of simplifying and learning from data, and also when
the dimensionality of the space is high.

To conclude, a great deal of space was dedicated to introducing and applying
PCA, as well as diagnosing correlations across features in the original high
dimensional space, because nearly every dimension reduction that has followed
has built on PCA. Regardless of the scope of any post-PCA methodological
innovation, PCA’s approach to making high dimensional spaces more manageable
and interpretable has not only impacting decades of development of subsequent
dimension reduction techniques but has contributed to the rise of an entire
subfield dedicated to dimension reduction. This subfield and its
interpretational value in high dimensional contexts has never been more
important than in the current era of massive data production at an incredible
rate.

###  2.5 Suggestions for Further Reading

Friedman, Hastie and Tibshirani (2001) present excellent detail on PCA
estimation in Chapter 14.5, along with discussion of extensions. Zou, Hastie
and Tibshirani (2006) offer an extension of PCA in high sparsity contexts,
whereas Schölkopf, Smola and Müller (1997) derive a nonlinear version of PCA.

##  3 Locally Linear Embedding

PCA is very valuable as a starting place for reducing complexity of some
higher dimensional data space and learning an orthogonal representation of
that data space based on maximizing variance across the input features. The
result is a clearing of the clutter introduced by multicollinearity across the
features. Yet, as Goodfellow et al. (2016) point out, we must move beyond a
simple linear combination of features if we want to more efficiently handle
complex, high dimensional spaces. In short, we need PCA to understand
dimension reduction, while launching to construct more complex algorithms to
help us handle more complex data spaces.

Locally linear embedding (LLE), which has been around for about 20 years, is a
linear dimension reduction technique like PCA (Roweis and Saul, 2000). Yet,
LLE’s focus on local structure and then global projection allows it to handle
and accurately summarize more complex and nonlinear data spaces. As PCA is
interested in fitting multiple summary lines to data in orthogonal directions
of greatest variance, LLE is interested in learning the shape of the
underlying data structure. This shape is called a manifold, and LLE proceeds
by linearly combining weighted features (like in PCA), but in a fundamentally
different way as we will soon see.

As such, I introduce several new concepts and themes to deepen coverage of
dimension reduction: nonlinear dimension reduction, manifold learning, local
vs. global approximation and representation. These naturally flow from the PCA
approach, and also set up the more complex approaches covered in the following
sections (e.g., we have to know what a manifold is if we want to use LLE,
UMAP, or t-SNE).

###  3.1 Manifolds and Complex Structure

To effectively reduce the dimensionality of data, we need algorithms that,
first, appreciate this reality of underlying data structure, and then, are
flexible enough to capture it and then project it locally in a simpler
subspace as we have noted to this point. A different way to think about data
complexity is in terms of a manifold. A manifold is a $d$-dimensional
geometric shape that is treated as locally Euclidean. Assuming a manifold
underlies data, then, means that regardless of the contours of the manifold
(that is, whether or not it is highly nonlinear), all observations are assumed
to be exist somewhere along it, and can often be treated locally as linear.
This unlocks the potential of applying many methods to learn the shape and
contours of the manifold in a higher dimension, which can then be approximated
and projected onto a simpler subspace. Though our goal remains the same, the
key difference by introducing manifolds is that we are no longer interested in
directions of variation across the features (i.e., data summaries), but rather
we are interested in learning the full nuance of the structure along which the
data are distributed, and then recreating that learned structure in a lower
dimensional setting.

Manifold learning is a central task in machine learning, and especially
unsupervised learning. The notion of a latent manifold has been around for
well over a hundred years, e.g., Riemann (1873). Indeed, Riemann’s ideas have
influenced not only the field of differential geometry but have given rise to
a subfield known as “Riemannian geometry.” At a basic level, Riemannian
manifold learning is premised on the idea that a latent, simpler manifold
characterizes the input space, but in some high dimensional ambient space.
Ambient in this context means the possible space an observation could occupy
is massive. Yet, that observation only occupies a single point in space. For
example, members of Congress represent districts, cities, and states. So, at a
strictly theoretical level, a legislator could be from the 35th district
representing San Francisco, Louisiana. There is a district number, a city, and
a state. The combinations are non-sensical in practice, of course, but this
combination is technically possible in ambient space. And thus, the ambient
space is comprised of a massive amount of possible places to occupy in the
high dimensional setting (35 in our running ANES example, where this silly
example is only a 3-dimensional case). Yet, as each observation only occupies
a single position across all dimensions in the higher dimensional, ambient
space, the assumption is that all input data are distributed in this space and
can be represented in a lower-dimensional way. The task of dimension reduction
in these terms, then, is comprised of two parts: first, learn the contours and
shape of the latent manifold (including distances between observations that
exist in the common space), and then second, recreate the true structure in a
lower dimensional, understandable way. In so doing, for these types of
dimension reduction tasks, we are moving from some massive high dimensional
ambient space with nearly infinite configurations across all dimensions, to a
still-high dimensional, but more real version of the data (full, raw input
space, $\mathbf{X}$), to ultimately a low dimensional subspace of two or
possibly three dimensions.

We can make this come alive a bit by demonstrating LLE’s ability to unravel a
three-dimensional “S-curve” of data, and project it to a two-dimensional
scatterplot. We will use colors in the following illustration to help with
understanding the location of points in the higher and lower dimensional
versions. Of note, this application is possible and requires little code, when
loading a few packages that contain the S-curve data. Start with the three-
dimensional S-curve data in Figure 3.2. To learn this structure and project it
onto a lower dimensional subspace, we would be interested in a two-dimensional
scatterplot of the data, but with the original structure remaining intact.
This projection is shown in Figure 3.2.

Figure 3.1: 3D S-Curve

Figure 3.2: 2D S-Curve Projection

In Figure 3.2 we can see the relatively complex S-curved structure. If we
wanted to flatten this curve so to speak, and plot it in two dimensions on a
scatterplot, a way to do so would be to first learn the distances between each
point across all three dimensions. Recording these dimensions as new
coordinates, we can then plot these coordinates in a two-dimensional sub (new)
space. This is what LLE is doing at a high level, and indeed, LLE was used to
flatten the S-curve in Figure 3.2 and then project the lower dimensional
result in Figure 3.2.

Critical to understanding dimension reduction from a manifold setting
perspective is the balance between local and global structure. Local structure
refers to the pointwise approximations of a space. Combined across all
observations in the data space, we then get a global picture of the true
structure of the data, based on the local, pointwise approximations. But
importantly, as with a photograph, the projection is not a perfect mirror of
the original high dimensional input space. There are errors that are made in
attempting to take a 35-dimensional set of coordinates and project it down to
a two-dimensional set. Yet, as with PCA, which was also based on “making
mistakes” by throwing data away for the benefit of focusing on the majority of
the variance across the features, we are willingly “throwing data away” in an
implementation of LLE in an effort to learn the contours of the latent
manifold. Once learned, we are then able to represent it in a way that is more
intuitive and understandable.

###  3.2 Formalizing LLE

LLE attempts to reconstruct a data space based on a series of $n$ pointwise
local calculations of similarity. Local means proceeding across the input
space on a point-by-point basis, and calculating an approximation of each
observation $i\in\\{1,\dots,n\\}$, on the basis of a weighted, linear
combination of a set of nearest neighbor observations. This approach to
locally learning the data space as it is distributed along the manifold allows
for the full space to be reconstructed, which results in the lower dimensional
projection. Importantly, LLE is a linear approach to nonlinear dimension
reduction. That is, the process of deriving the weighted local representations
of each point, which are later combined to represent the global structure, are
linearly combined (i.e., added together). Yet, the shape and structure of the
manifold along which the data are distributed can be (and often is) nonlinear
and complex, especially as the dimensionality of the space increases. The goal
of LLE, then, is to find the weights in a local fashion, combine (multiply)
them by the raw feature values, and then add them together to give an
approximation of the single point. This process is repeated for all
observations, up to $n$ times, resulting in the learned, lower dimensional
representation of the high dimensional input space.

To formalize LLE, start with the high dimensional input space, $\mathbf{X}$,
which recall is an $N\times P$ data matrix, consisting of elements
$i\in\\{1,\dots,n\\}$ across column features, $p\in\\{1,\dots,P\\}$. The
output space, which given our focus in this Element, is not a predicted value
like in a supervised task, but is another matrix, which we can call
$\mathbf{Y}$ for simplicity. We let $\mathbf{Y}$ be an $n\times d$ matrix,
consisting of values $i\in\\{1,\dots,n\\}$ across each dimension,
$d\in\\{1,\dots,D\\}$. Importantly, as with PCA and the other techniques
covered later, we will typically let $D=2$, such that we are interesting in
moving from the higher dimensional space, e.g., in our case $D=35$, to a lower
dimensional subspace, e.g., $D=2$. This is a formal way of describing any two
dimensional plot like the earlier PCA scores plot in Figure 2.13.

Recall we are interested in linearly combining weighted with original feature
values to produce an approximation of each point, ultimately done $n$ times.
These weights, which are at the heart of LLE and will be stored in
$\boldsymbol{\beta}$, are based on distances between a candidate observation,
$i$, and a small surrounding neighborhood of observations, $j$, of size, $k$.
Thus, the weights are subject to a normalizing factor for each observation,
$i$,

$\sum_{j=1}^{k}\boldsymbol{\beta}_{ij}=1.$ (3.1)

Regarding distance, most applications of LLE and other dimension reduction
algorithms either use Euclidean distance or Manhattan (“city block”) distance.
Euclidean distance, $d_{e}$, is found by taking the square root of the squared
difference between observations $i$ and $j\forall i\neq j$. For example, in a
single dimension, which is easily scalable to higher dimensions,

$d_{e}(i,j)=\sqrt{(i-j)^{2}}.$ (3.2)

Manhattan distance, $d_{m}$, measures the distance between two points by
taking the absolute value of the difference between two observations, $i$ and
$j\forall i\neq j$. Again, in one dimension we have,

$d_{m}(i,j)=|i-j|.$ (3.3)

For our purposes, unless otherwise noted, whenever a measure of distance is
required in the techniques covered in the Element, Euclidean distance will be
used.

Imagine a small neighborhood of observations surrounding each candidate
observation, $i$. We are interested in calculating the distance between $i$
and each surrounding observation, $j$. We then multiply the distance, which
itself is the “weight,” by the raw feature value for each feature across all
points in the small neighborhood. Importantly, allowed by the constraint in
Equation 3.1, we can more readily interpret the weights for each point, $j$,
relative to $i$, such that higher values indicate the observation is closer to
and thus looks more like the candidate observations, $i$. In this way, across
$n$ points, the local structure is learned in a point-by-point fashion and is
then able to be represented in a lower dimensional version of the space, as we
are left with a set of $n\times d$ coordinates (stored in $\mathbf{Y}$) by
which we can simply produce the $d$-dimensional plot. LLE accomplishes this by
capturing and retaining all of the original information and structure across
the full data space (which again, is distributed along the latent manifold).
By plotting $\mathbf{Y}$, we get the simpler, lower dimensional representation
of the input space that has retained the original, high dimensional structure.
This is how LLE can produce Figure 3.2 given only the complex input data in
Figure 3.2.

Of note, LLE, as with many of the other techniques we will cover like
autoencoders, is all about reconstruction. That is, there are many possible
combinations of neighbors that could surround $i$, and be used to approximate
that point. Yet, the optimal LLE solution is defined by the configuration of
all local approximations that minimize the reconstruction error. As
reconstruction is based on the approximated value of $i$, based on the
weighted sum of the neighborhood of $j$ of size $k$, this suggests we need to
minimize the mistakes we make in making the pointwise calculations. Put
differently, reconstruction error in the context of LLE is identical in
structure to error in a simple linear regression setting. As in regression and
LLE, our goal is to find the solution that minimizes error. Error in both LLE
and regression settings is defined as information loss. By approximating the
relationship between the inputs, $\mathbf{X}$ and the output, $\mathbf{Y}$, we
are necessarily losing information for the sake of a well fitting, but
parsimonious explanation of the data. To demonstrate this point, consider the
similarity in the loss functions for regression in Equations 3.4 compared to
LLE in Equation 3.5.

For regression problems, the most common approach to define loss is the mean
squared error,

$MSE=\frac{1}{N}\sum_{i=1}^{N}{\left(y_{i}-\hat{f}(x_{i})\right)^{2}}.$ (3.4)

The goal in regression is to identify a model that generates the smallest
possible MSE. That is, we want to minimize the error (mistakes in prediction)
when we decide to approximate some relationship between in the inputs and
output.

The goal is the same in LLE. By substituting $y_{i}$ in Equation 3.4 with the
candidate observation, $x_{i}$, we are treating the candidate observation as
the “output” for a single step in the algorithm. With a slight update to the
notation, we get a measure of reconstruction error for a fit of LLE,

$RSS=\frac{1}{N}\sum_{i=1}^{N}{\left(x_{i}-\sum_{\forall j\in k\neq
i}\beta_{ij}x_{j}\right)^{2}}.$ (3.5)

By using the set of weights found from minimizing RSS in Equation 3.5, we are
able to find the set of two-dimensional coordinates, which ensures a minimal
global reconstruction error, based on use of the optimal set of weights.

The optimal LLE solution, then, is found when RSS is smallest. By locally
approximating each point in such a way, we can find the set of new
observations, $\mathbf{Y}$, which are derived on the basis of the weighted
contributions of each feature, $p$, each point $j\neq i$, across all
observations, $n$, distributed along the manifold. The result, when plotted,
should look similar to all original candidate observations around which each
neighboring observation was based, only now in two dimensions rather than 35.
As such, LLE is able to capture global structure on the basis of local,
pointwise approximation.

###  3.3 Applying LLE to the ANES Data

As with the section on PCA, we will transition to apply LLE to the ANES 2019
pilot study data. We will be using the same input features covering feeling
thermometers on a battery of issues, people, countries, and institutions.
Prior to exploring the data, we first load the appropriate packages and scale
the data to ensure all data are effectively unitless, and thus comparable on a
common scale.

With the data loaded, we will explore the data in a similar 3D way as in
Figure 3.2. We will explore three features at a time that are substantively
related to each other to begin to get a sense (albeit an incomplete one) of
the shape of the manifold. This is necessarily “incomplete” because a
35-dimensional plot would be substantively useless. See the results in Figure
3.3.

Figure 3.3: 3D Distributions Across 12 Features

The four plots in Figure 3.3 show the distributions of three features unique
to four substantive topics: Politicians (Bernie Sanders, Donald Trump, Barack
Obama), Institutions (the UN, NATO, the NRA), the issue of immigration (ICE,
illegal immigrants, immigrants in general), and middle eastern countries
(Israel, Palestine, and Saudi Arabia). These plots in Figure 3.3 reveal some
useful patterns to inform our LLE application. For example, respondents tend
to have extreme feelings toward the politicians with clear groupings of
respondents at the bounds of the feeling thermometers on all three axes.

We need to search for the optimal neighborhood size, $k$, prior to fitting the
model. The main lle package includes a useful function, calc_k(). This
function runs LLE for each supplied value of $k$, which in the case below was
a real-valued, non-negative integer from 1 to 20. The “optimal” value is
defined at the value of $k$ for which the coefficient, $\rho$, is highest.
$\rho$ is simply the correlation between approximations in the high
dimensional setting compared to the low dimensional projection. Higher
correlation means better representation. The function calc_k() allows for
plotting values of $1-\rho^{2}$, such that we are interested in the value of
$k$ for which this quantity is lowest.

Using parallel computing on two different computers, it took 9.2 minutes run
on 7 cores, and 10.9 minutes to run on 3 cores. Within the function call, we
first pass the scaled data object (minus “party affiliation”), followed by
specifying the number of dimensions for the projection, $m=2$. I then set
parallel = TRUE and specified the number of cpus to use for the parallel
process. Parallel programming of this sort is not necessary, but it
significantly speeds up the code, especially for larger data applications.

Figure 3.4: Optimal $k$ from Minimal Value for $1-\rho^{2}$

In addition to the code producing Figure 3.4, I include a brief line to
manually locate the optimal value of $k$. Both the figure and the manual
search reveal the optimal neighborhood size to give the best representation of
the high dimensional space is 19. We proceed accordingly, setting $k=19$. We
follow the run of the algorithm with a visualization of the coordinates that
mirrors the PCA approach. See the results in Figure 3.5.

Figure 3.5: LLE Projection Results

From Figure 3.5, there indeed seems to be natural partisan differences across
the first two dimensions in the projection space. Recall, the feature for
party affiliation was not included in the LLE fit, allowing us to
conditionally color points accordingly. Though differences on a party
dimension seem to exist, the LLE results can be better contextualized by
directly comparing with scores from a PCA fit. To do so and directly compare
patterns, see Figure 3.6.

Figure 3.6: Comparing LLE to PCA

First, in Figure 3.6, we can see similar distinction between Democrats (blue)
and non-Democrats (red). Both algorithms are picking up on the partisan split
in the feeling thermometer space. Interestingly, though, PCA results in a more
diffuse projection, whereas the LLE results place respondents in a more
compact space. As the LLE algorithm is capturing local structure and
projecting the nuance of the distribution of points along the manifold,
compared to PCA, which is plotting the most unique variance in the data along
the first two dimensions, the LLE results are likely more accurately picking
up the shape of the underlying structure of the data. Whereas the PCA solution
is helpful to show the distribution of variance across the input space in the
first two dimensions, the LLE solution is likely giving a more realistic look
at the overall shape of the structure. This suggests, then, that the true
distance between Democrats and Republicans, though distinctly separate, is not
as diffuse as the PCA solution seems to suggest. The contours of this
separation are addressed throughout the coverage of the other techniques in
the Element.

To recap, the goal of LLE is to create a simpler representation of the complex
original space based on reconstruction error (and an attempt to minimize it).
As such, the idea of reducing the original space with an eye toward recreating
it on the basis of the simplified version is central to another approach to
nonlinear dimension reduction called autoencoders. Autoencoders, which are
built on a neural network-based architecture, are addressed later in this
Element.

Importantly, LLE helps address the limitations of PCA by focusing on recovery
of the latent structure through both local and global exploration of the data,
rather than on simply maximizing shared variance as in PCA. Yet, a weakness of
LLE is the linear combination of the weighted features, which can limit LLE’s
ability to recover global structure, especially in non-convex contexts for
example. While LLE can handle nonlinear data better than PCA, it cannot do so
as efficiently as some of the more recent, graph-based (t-SNE or UMAP) or
neural network-based (SOM or autoencoders) techniques covered in the following
sections.

###  3.4 A Risk-Averse Workflow

To this point, we have covered PCA and LLE, along with several initial checks
for multicollinearity across the feature space. As such, this final subsection
offers a very biased opinion on what a risk-averse workflow for a dimension
reduction project might look like given weakness of and connections between
PCA and LLE.

First, adopt a demeanor of extreme caution when pursuing and interpreting
output from dimension reduction algorithms and models. As these are
unsupervised, it is impossible to conclude that the patterns uncovered are
indeed a good and reliable reflection of the high dimensional concept you are
trying to capture. With that, though, I recommend starting with PCA as it is
widely used, statistically well-established, and interpretation is quite
intuitive. Upon generating the PCA solution, plot the scores as these are
essentially the new measurements for the newly calculated features/components.
These are the features that would be fed to a supervised or other model
downstream, depending on the goals of the project. And to help deepen
practical understanding of the PCA solution, consider coloring the scores in
the two-dimensional plot along some theoretically interesting feature as in
Figure 2.13. Then, with the single plot in hand, proceed to a more complex
algorithm like LLE and compare the patterns as in Figure 3.6.

A final step that never hurts is to compare to the original, raw feature
values. That is, plot several input features against each other and still
color by party affiliation (or the key conditional feature) to explore whether
natural separation in the high dimensional space is similar to the recovered,
low dimensional separation. Importantly, as this is a manual, feature-by-
feature approach, it is inefficient to progress across all input features, and
indeed may not make sense to do for especially complex, high dimensional data
sets (e.g., $p>100$). Yet, for the sake of demonstration, we might plot a few
key features against each other and color densities of observations by party
to explore whether we find similar separation in the original space as we do
in the transformed, reduced space as in the plots in Figure 3.6. To do so, and
merely for demonstrative purposes, I plot the two-dimensional densities of
observations across eight features, and color by party affiliation. They are
feelings toward: Trump and Obama (upper left in Figure 3.7), ICE and illegal
immigrants (upper right in Figure 3.7), the UN and NATO (lower left in Figure
3.7), and Palestine and Israel (lower right in Figure 3.7). See all of these
results in Figure 3.7.

Figure 3.7: 2D Contour Plots Across Eight Original Features

In Figure 3.7, we can certainly see separation in the raw, original feature
space. For example, in the Obama/Trump plot in the upper left of Figure 3.7,
observations of each party are extremely far from each other as expected
(i.e., feelings toward Trump are extremely polarized relative to feelings
toward Obama). We see similar patterns for feelings toward ICE and illegal
immigrants, and to a lesser degree in the feelings toward Palestine and
Israel. Yet, perhaps most interesting is the spread of raw feelings toward
NATO and the UN in the lower left plot of Figure 3.7, which reveals a
strikingly similar pattern of separation that we saw in the projection space
from both PCA and LLE in Figure 3.6. There is a cluster of overlapping
respondents in the middle of the plot, but large and clear clusters at the
outer bounds across respondents of different parties. As a result, perhaps
these features could be proxies for broad partisan divisions between
respondents.

This workflow based on comparison from several angles (e.g., across algorithms
and across original raw features) helps with interpretation. Patterns from
PCA, which recall are more interpretable than more complex manifold-based
learners, can be directly compared and thus placed into context. If patterns
from the more complex algorithm look widely different from the patterns in the
simpler PCA approach, then this could be a signal that something is less clear
than anticipated. And further, it could provide a useful clue to you, the
researcher, that more careful thought is required to effectively learn the
underlying structure of the data. If, on the other hand, similar patterns
emerge, then this lends support for your expectations that, first, some non-
random structure does indeed appear to characterize these data, and second,
you have homed in on its contours and shape. That is, you have deepened an
understanding of your data as you set out to do in this unsupervised,
exploratory type work. Indeed, this is the general approach I take in this
Element. I will continually compare substantive patterns with conditional
coloring in my plots as a sanity check, of sorts, to ensure the patterns I am
uncovering are likely non-random, and real.

Ultimately, though, as referenced frequently to this point and elsewhere in
this Element, we are engaging in unsupervised machine learning, such that no
ground truth is guiding the process. As a result, we are forced to proceed
with our best judgement, which is by definition always ambiguous to at least
some degree. The best we can do is compare across several approaches to a
common problem using common data, and critically assess the emergent patterns
to give insight into our data and motivating problem.

###  3.5 Suggestions for Further Reading

Building on the original LLE derivation (Roweis and Saul, 2000), Saul and
Roweis (2003) present an explicit machine learning view of LLE, which is
useful to understand more complex extensions including the robust version of
LLE offered in Chang and Yeung (2006).

##  4 Nonlinear Dimension Reduction for Visualization

As we have seen to this point, we can use dimension reduction techniques to
reduce the complexity of data, giving a simpler, more manageable version of
the high dimensional input space. This is especially useful when the size and
dimensionality of the data are quite large, which is less common in social
science applications. Even still, as we have begun to see, the precise
approach to reducing the dimensionality of a high dimensional space, can vary
quite widely. We can consider a simpler version of a data space as one that
reflects the most variance across the input features (PCA). Or we could
consider a simpler version of a data space as one that reflect local and
structural similarity, giving a lower dimension approximation or the original
space (LLE).

Another way to think about how to select from among the many dimension
reduction techniques is based in the goal of the task or project. To this
point, we have been implicitly interested in reducing the complexity of a data
space to result in some reduced set of newly-calculated features. Though this
process of feature extraction, we proceed to plot these new features against
each other and then condition the color of the observations by their party
affiliation to offer an informal look at validating the quality of the
dimension reduction solution. That is, did the technique separate members of
different parties on the basis of a simplified version of the feeling
thermometer space?

In line with the approach to this point of plotting results, is to more
explicitly treat dimension reduction as a tool for visualization. That is,
reducing the complexity of a high dimensional data space with an eye toward
assessing the patterns visually, rather then an a priori focus on the new sets
of calculations giving the new features as in PCA and LLE. Two extremely
popular and widely used techniques for this goal of dimension reduction for
visualization are: t-distributed stochastic neighbor embedding (t-SNE) and
uniform manifold approximation and projection (UMAP).444Importantly, there
exist other approaches to dimension reduction for visualization such as
multidimensional scaling and ISOMAP. Though valuable, these are not covered in
this Element due to limitations in space. For excellent coverage of
multidimensional scaling in the social sciences, see Armstrong et al. (2014).
In this section, we will cover both t-SNE and the more recent UMAP approaches
to dimension reduction focusing especially on the visual results of the
algorithms to help us move closer to out goal of understanding whether latent
partisan trends exist in the higher dimensional feeling thermometer space.

###  4.1 From Linear to Nonlinear Dimension Reduction

Though the task focused on in this section is only slightly distinct from the
previous sections and though all approaches have included heavy visual
components, a key shift in covering t-SNE and UMAP is in how the algorithm
treats the data to find the optimal lower dimensional projection. Recall, PCA
and LLE were both linear approaches to dimension reduction, finding the lower
dimensional data space based on adding together a linear combinations of
weighted raw feature values. With t-SNE and UMAP, though, we rely on a
nonlinear combination of raw input features to give a very different treatment
of the data. We are still working with the same high dimensional data (35 to
be precise), but we are now treating and processing it a bit differently than
we have to this point.

Both t-SNE and UMAP rely on neighbor-based smoothing to give local versions of
the high dimensional data. That is, we are interested no in working locally
and in a neighbor-based fashion (as with LLE), but based only on distances
within the neighborhoods, rather than using those distances to give a weighted
version of the raw input features. Then, the manifold is attempted to be
recreated based on these local distances. The distances are calculated such
that observations far away are also pushed far away in the projected version
of the full space, whereas observations that are close to the candidate
observations in the high dimensional setting are brought very close to the
candidate observation in the lower dimensional project. These repulsive and
attractive forces are baked into the algorithm, which result in a lower
dimensional version of the high dimensional space with spatial similarities
and differences being exaggerated. Spatially similar observations are grouped
tightly together, whereas spatially different observations in the higher
dimensional setting are reproduced to be pushed far from each other.
Importantly, t-SNE is a probabilistic approach to this problem whereas UMAP is
a graph-based approach that relies only on smoothed distances metrics. This
innovation in UMAP overcomes the core limitations in t-SNE, which is addressed
in the remainder of this section.

The result from these approaches to dimension reduction is a (usually) two-
dimensional plot of the lower dimensional version of the data space that
looks, at first glance, quite odd given these exaggerated differences between
observations on the basis of spatial similarity and difference. Yet, with the
ability to tune the key hyperparameters in each algorithm to control the
global versus local behavior of the algorithm (i.e., the amount of
exaggeration), these algorithms can give incredibly important insight into the
underlying structure of high dimensional data.

###  4.2 t-SNE

t-SNE is an algorithm that was developed in 2008 in an effort to make clear
the underlying structure in data in a visually digestible way. With an
expressed focus on developing and using the algorithm for visualization,
Maaten and Hinton (2008) center their approach on comparing probability
distributions. The comparison between the original high dimensional version of
the data is compared to the new lower dimensional representation using the
common distribution-comparison metric, Kullback-Leibler (KL) divergence. When
values are small, this means the distributions are extremely similar, and when
KL-divergence is high, then distributions are different from each other. But
first, what are these distributions being compared?

The first step of t-SNE is to measure the distances between each candidate
observation, $i$, and the surrounding neighbors, $j$ in the neighborhood the
size of $k$. These distances are then placed under a normal distribution in
order to convert them to probabilities, where higher probabilities mean an
observation, $j$, is close to observation $i$. To make these more
interpretable, $d(\cdot)\forall j\in k$ for each $i$ are scaled to sum to $1$
and can be stored in $\mathbf{Z}$, where $d(\cdot)$ is a measure of distance
between some $j$ and the candidate $i$,

$\sum_{j=1\forall j\neq i}^{k}\mathbf{Z}_{ij}=1.$ (4.1)

The set of coordinates capturing the probabilities between each observation in
the high dimensional setting is stored in $x_{h}$, which belongs to the full
probability space $\mathcal{X}$.

Then, the next step is to create a new, randomly-distributed low dimensional
version of the data set in (usually) two dimensions. That is, the algorithm
randomly plots the data in a lower dimensional setting, where $d=2$. Then, the
process of placing a distribution around each observation, $i$, and the
calculating probabilities of observations lying close to each other is
repeated. Yet, this time, the t-distribution is used to account for greater
uncertainty given the force of information loss by placing the similarities
into a lower dimensional setting, i.e., moving from $d=35\rightarrow d=2$.
With a new set of probabilities between all observations and each other in
hand, the moving of observations continues at each iteration of the algorithm,
continuing to move observations closer together that have increasingly
probability of being close to each other, and also moving farther apart from
observations that have a lower probability of being similar to each other. We
are left with a lower dimensional version of the data space
$x_{l}\in\mathcal{X}$.

Finally, we come to the distribution comparison. At each iteration (the
movement of data points to place similar observations closer and dissimilar
observations farther apart), the algorithm compares the matrix of
probabilities between all observations in the lower dimensional space,
$x_{l}$, to the original matrix of probabilities in the original, higher
dimensional space, $x_{h}$. The measure of cost that captures differences
between these versions is the KL-divergence metric, which is sometimes called
relative cross-entropy,

$\sum_{x\in\mathcal{X}}x_{h}\text{log}(\frac{x_{h}}{x_{l}}),$ (4.2)

where each probability distribution $x_{*}$ in the full probability data
space, $\mathcal{X}$, is compared by logged, relative difference between the
high dimensional version, $x_{h}$ and the low dimensional version, $x_{l}$.

The idea here is that as probabilities in the lower dimensional setting look
more like probabilities across observations in the high dimensional setting,
then the KL metric decreases, suggesting less dissimilarity in the data space.
When the metric is low and unchanging, the algorithm stops, suggesting we have
arrived at a lower dimensional version of the high dimensional data, that has
a structure that mirrors the original complex structure. And to reiterate, we
are no longer working with weighted versions of the original input features.
Rather, we are now attempting to capture the probabilities that observations
lie close to some in space and far from others in the same space. Basing these
measures of similarity and differences on measures of spatial similarity,
instead of linearly adding individual feature contributions together, we are
bypassing the additive component of LLE and PCA’s approach to dimension
reduction. This is what is meant by “nonlinear dimension reduction” in the
case of t-SNE.

Thus, the goal of t-SNE is to reproduce the probabilistic (via spatial)
similarities across all points in the high dimensional context (again, based
on local probabilistic similarities), in a lower dimensional (usually 2)
setting, such that clearly similar groups of observations are tightly
connected and are distant from clearly different groups of observations.

There are multiple hyperparameters that need to be tuned, and thus globally
applied to the model to give an optimal solution. The most important
hyperparameter is perplexity. Perplexity controls the tradeoff between focus
on reproducing global structure versus local structure. Ideally the goal is to
effectively balance between these, thereby finding a good lower dimensional
representation of the more complex higher dimensional input space. Perplexity
controls this tradeoff by varying the size of the distribution placed around
each point. The higher the perplexity value, and thus the wider the
distribution around each point, the more global the solution will look, as the
distribution will allow point potentially far off in space to be given a
higher similarity/probability of being close to the center of the
distribution, which is where the candidate observation, $i$ is. So, by
increasing the size of perplexity and thus the size of the pointwise
distributions, we are allowing for higher a greater probability that distant
points are treated as more similar. Inversely, then, lower perplexity gives
much narrower distributions around each point, resulting in a lower
probability that distant points will be treated as similar to the candidate
observation, $i$, at the center of the distribution. Thus, lower perplexity
scores result in recovery of more local structure. The question becomes, how
do we home in on the optimal perplexity value, given the potentially drastic
change in the lower dimensional representation we get as a result? There are
two options here: manually vary the value and inspect distributions or conduct
a grid search across multiple parameter values to inspect the evolution of
structure over the range of the hyperparameter.

Consider the first approach of manually varying perplexity. For the
application of t-SNE using the ANES data, we will rely on the Rtsne package.
As before, we will be working with the ANES clean data object and several
tidyverse tools (e.g., ggplot2, dplyr). Upon manually varying perplexity, we
can patchwork the six plots together to result in Figure 4.1.

Figure 4.1: Manually Varying Perplexity for t-SNE

Note the clear change in the picture of structure we get in Figure 4.1, when
perplexity is varied manually. As expected, moving from small values of
perplexity we pick up on more local nuance in the structure, compared to high
perplexity values, where the global structure is the focus. In the case of the
latter, e.g., perplexity $=500$ as the extreme version in the lower right
plot, the space is extremely dense, suggesting “structural similarity,”
essentially defined by the fact that the respondents were simply responding to
the same survey. In other words, such extremely high perplexity yields a less-
then-useful look at this space, where we are unable to detect any degree of
separation in the projection space.

Regardless of the different values of perplexity, a clear pattern is the
separation between Democrats and Non-Democrats in this projection space.
Though the separation naturally lessens when we tune the algorithm to focus
more on global structure, the parties are never fully mixed together. This
suggests that there is indeed a latent partisan structure to this feeling
thermometer space.

The Rtsne package is built on a subtly updated version of t-SNE using the
Barnes-Hut approximation. Though details of this extension are beyond the
scope of current purposes, I point interested readers to the details in Van
Der Maaten (2014). In brief, this implementation of t-SNE includes an
additional hyperparameter, $\theta$ to capture the tradeoff between speed of
fitting of the algorithm with the accuracy of the solution. Therefore, to more
efficiently search across both of these hyperparameters, the code below
briefly demonstrates that which such a grid search might look. Importantly,
though limited by space, the following code can be extended in numerous ways,
e.g., comparing across other hyperparameters such as $\eta$, which captures
the learning rate or size of steps to take between iterations. But sticking
with perplexity and $\theta$, for which zero gives the original t-SNE
implementation and thus the most accurate, but slowest solution, whereas
$\theta=1$ gives the fastest, but less accurate solution. Of note, as we will
see in the follow subsection, UMAP is a much faster algorithm being run in a
fraction of the time as t-SNE. Indeed, the t-SNE grid search took about 19
minutes to run on 4 cores.

Figure 4.2: t-SNE Grid Search Across $\theta$ and Perplexity

From Figure 4.2, there is a clear shift in the structure across various
hyperparameter values. For example, as previously seen in Figure 4.1, higher
perplexity values reveal a more global version of the lower dimensional
projection, whereas smaller perplexity values in the left-hand plots of Figure
4.2 pick up on more local, diffuse structure. Yet, interestingly, not much
nuance appears for either lower or higher $\theta$ values, across the range of
perplexity. The clear shift in structure is seen in the columns as perplexity
is varied. The goal of this demonstration is offering a more efficient
approach to varying multiple hyperparameters in a single chunk of code,
instead of manually varying values and checking patterns after each.
Especially for a slow algorithm like t-SNE, the manual approach is much less
efficient.

###  4.3 UMAP

At this point, we come to a close relative of t-SNE: UMAP (McInnes, Healy and
Melville, 2018). UMAP is a very recent approach to nonlinear dimension