
\chapter{Data Visualization\label{chapF:Directions-for-Visualization}}

This chapter considers OODA \index{visualization|(}visualization
concepts in more detail. Here the focus is on data objects having
a Euclidean feature space\index{feature space} representation, where
the entire data set is conveniently summarized by a matrix. For the
reasons given in Section \ref{secBA.1:Foundations_Terminology}, columns
of the matrix are taken to represent data objects. A direct and natural
visualization of the structure of a matrix of data is a \emph{\index{heat-map|(}heat-map}
image, whose application in OODA contexts is studied in Section \ref{secF.A:Heatmap-Views}.
\emph{Curve} and \emph{scatterplot matrix\index{scatterplot matrix}
}views\index{scores!scatterplot|(} of data objects, both of which
have been used repeatedly in earlier chapters, are developed in further
detail in Sections \ref{secF.B:Curve-Views} and \ref{secF.C:Scatterplot-Matrix-Views}.
Data centering and combined views are considered in Section \ref{secF.E:ComboView-Centering}.

\section{Heat-Map Views of Data Matrices\label{secF.A:Heatmap-Views}}

The main idea of a heat-map view of the structure of a data matrix\index{data matrix!heatmap view}
is to construct an image that represents each matrix entry with a
colored (or gray level) patch in a rectangular grid. Patterns in the
perceived colors frequently give useful insights about the data set
represented by the matrix. The basic idea is related to that of digital
photography (where the patches are called \emph{pixel\index{pixel|(}s}
from picture elements), and has been used in many contexts. A good
overview of the area, that includes some interesting historical discussion,
can be found in \citet{Wilkinson2009Heatmap}. This approach to data
visualization is commonly used in bioinformatics, where it has proven
to be very effective in finding biological insights and discoveries
from data, see for example \citet{Eisen1998ClusterHeatmap} and \citet{perou2000molecular}
where this exploratory\index{exploratory analysis} tool was used
in the discovery of clinically relevant cancer subtypes, that have
since become key to realizing the improved treatments promised by
\emph{precision medicine\index{precision medicine}}. However, heat-map
visualization has also attracted skeptics such as \citet{Rogowitz1996NotLieViz}
and \citet{Borland2007RainbowColors}, particularly over color choices
and how they relate to human visual perception\index{visual perception}.
A different limitation of heat-map data visualizations is shown in
Figures \ref{figF.C1:ToyHeatMap} and \ref{figF.C2:ToyPCA}, where
it is hard to see small signals in noisy data. For visualization of
a matrix with nonnegative entries a gray level plot, as shown in Figure
\ref{figF.A:ToyHeatMapCluster} below, can be a good choice because
the ordering of colors is clear. A richer color scheme in the nonnegative
case is the topographical colors\index{topographical colors} that
are classically used in geographic maps: green, yellow, orange, brown,
gray and white. For entries that are both positive and negative, in
contexts where knowing which entries are 0 is important, two colors
plots such as blue and red are usually more useful, as discussed around
Figure \ref{figF.C1:ToyHeatMap} below. The topographic color scheme
can also be extended to this case by including shades of blue suggesting
below sea level. In situations where location of the origin is not
important, e.g. log scale data such as the Spanish Mortality data\index{Spanish Mortality data}
in Figures \ref{fig:Mortality-Raw}, \ref{figD-L:TCGA-Log2vsRaw}
and a number of other places in this monograph, a rainbow color scheme
can also be useful.

It is usually important to keep in mind that ordering of the rows
and columns in heat-map visualizations has a major influence on what
can be learned. This is illustrated in the \index{Two Clusters data}Two
Clusters data toy example in Figure \ref{figF.A:ToyHeatMapCluster}.
All 3 panels are based on a common set of values. At first sight the
left hand panel appears to be fairly random, with perhaps hard to
interpret vertical and horizontal patterns. In the center panel, a
clustering algorithm (\index{hierarchical clustering|(}hierarchical
clustering with Euclidean distance and average linkage\index{linkage!average|(}
as discussed in Section \ref{secK-B:HierarchicalClust}) has been
applied to the columns, and they are correspondingly reordered. Despite
the high noise level this shows a clear propensity for lighter shades
of gray in the right half. The right hand panel shows the results
of similarly clustering and reordering the rows as well, which shows
a related pattern between the top and right half of the plot. The
contrast between the left and right sides indicates two clusters in
the data set (columns) with different cluster mean structure evident
in the reordered variables (rows). This ability to show rich structure
(even in the presence of very high noise) seems to be why this type
of view is commonly used in bioinformatic contexts, where there tend
to be strong underlying patterns often corrupted by a relatively high
noise level. 

\begin{figure}[H]
\includegraphics[scale=0.9]{ChapterF/OODAbookChpFFigAToyHeatMapClusterCMYK}

\caption{Toy example demonstrating usefulness of hierarchical clustering rows
and columns in heat-map image views. Left panel shows a random ordering.
Column clustering has been applied in the center panel, and both column
and row clustering were used in the right panel, to reveal important
systematic ``Two Clusters'' structure.\label{figF.A:ToyHeatMapCluster}}
\end{figure}
In some situations, the hierarchical clustering itself is of interest,
beyond merely enhancing the visualization. In such cases it is common
to add dendrograms (as shown in Figures \ref{figK-E:DendorEucSing}
and \ref{figK.F:HierHiDDendros}), revealing the tree of \index{nested analyses}nested
clusters for both the rows and columns of the heat-map.

Another issue that is critical to heat-map visualizations is \emph{scaling}
as illustrated in Figures \ref{figF.B1:ToyHeatMapScalingLinear},
\ref{figF.B2:ToyHeatMapScalingLog} and \ref{figF.B3:ToyHeatMapScalingQuantile}.
These are three different views of the same data matrix demonstrating
the impact of different gray scale choices on which aspects of the
data are highlighted. Figure \ref{figF.B1:ToyHeatMapScalingLinear}
shows how the typical first choice of an equally spaced color scale
can obscure important structure in the data. In particular, the only
thing visible is one small white spot in the heat-map shown in the
left panel. Insight into the cause of this behavior comes from studying
the distribution of matrix values, which is shown in the right panel.
This uses the same format employed for displaying one dimensional
distributions in many other contexts above, starting with the Spanish
Mortality data\index{Spanish Mortality data} in Figure \ref{fig:Mortality-PC1}.
In particular, the matrix entries lying in the interval $[0,20]$
are shown as ($2500=50\times50$) circles. The curve is a smooth histogram,
i.e. kernel density estimate\index{kernel density estimate} as discussed
in Chapter \ref{chapN:Smoothing}. The horizontal coordinate of each
circle is the entry of the matrix and the vertical coordinate is the
order in the \index{vectorize}vectorized version of the data where
the matrix columns have been successively concatenated into a single
long column. Both the circles and the density estimate reveal a highly
skewed distribution. The right hand panel also contains vertical dotted
lines which are the boundaries of the 20 equally spaced gray scale
regions used in the heat-map in the left panel. The right hand three
quarters of the plot has very few circles, which reflects very few
white or gray pixels visible in the heat-map. But note that the shape
of the distribution of the circles suggests there may be additional
structure that is not visible in the heat-map because all those points
appear in the essentially black region. While equal spacing of color
scales is a common default in heat-map software packages, this example
demonstrates how that choice may obscure important structure.

\begin{figure}[H]
\includegraphics[scale=0.9]{ChapterF/OODAbookChpFFigB1ToyHeatMapLinear}

\caption{Toy example illustrating gray scale color scaling issues. Heat-map
image is shown in the left panel. Distribution of the matrix entries
is shown in the right panel, with vertical dotted lines indicating
the boundaries of the gray color regions. The heavily skewed distribution
indicates poor use of the equally spaced gray levels, resulting in
an ineffective heat-map view. \label{figF.B1:ToyHeatMapScalingLinear}}
\end{figure}
A simple approach to making better use of the available gray scale
is to transform the distribution of matrix entries to have less skewness\index{skewness}.
The shape of the distribution shown on the right of Figure \ref{figF.B1:ToyHeatMapScalingLinear}
suggests a log transform could be useful. The effect of this, for
the data matrix in Figure \ref{figF.B1:ToyHeatMapScalingLinear},
is shown in Figure \ref{figF.B2:ToyHeatMapScalingLog}. Note that
the heat-map image on the left still shows the same white spot as
in Figure \ref{figF.B1:ToyHeatMapScalingLinear}, but now it also
shows three distinct gray spots. The distribution of the log matrix
entries in the right panel clearly shows the impact of the log transformation\index{transformation!log}.
The bright spot in the image on the left is reflected by the larger,
sparser peak in the matrix values. The less bright but light gray
spots in the image appear as the peaks extending about halfway across
the range (thus using only about half of the gray levels). There are
two such peaks because of the column reorganization operation used
to order the matrix pixel values. A careful look at the heat-map reveals
additional interesting structure. In particular a regular grid of
quite dark gray (versus the black background) points is visible. Because
these small dark gray spots appear all across the heat-map, their
pixels form a thick band on the left side of the right hand panel,
itself indicating interesting structure in this data matrix. 

\begin{figure}[H]
\includegraphics[scale=0.9]{ChapterF/OODAbookChpFFigB2ToyHeatMapLog}

\caption{Different heat-map image view of the same data matrix from Figure
\ref{figF.B1:ToyHeatMapScalingLinear} using gray level scale (boundaries
again shown using vertical lines in the right panel) based on log
of function values. Distribution of log function values shown in the
right panel, with vertical lines now showing better use of gray levels,
so this color scale reveals more structure in the image. \label{figF.B2:ToyHeatMapScalingLog}}
 
\end{figure}
While this log scale\index{scale!log}d heat-map image reveals much
more structure in the data matrix than was available from the naive
scaling shown in Figure \ref{figF.B1:ToyHeatMapScalingLinear}, there
is still room for improvement. For example the entire brighter half
of the range (thus the bright half of the gray levels) is devoted
to just the single tallest peak in the matrix. An alternate approach,
that can be very revealing in some contexts is \emph{quantile\index{quantile!scaling}
scaling} (also known as \emph{histogram equalization}\index{histogram equalization}).
That color scaling ensures approximately equal numbers of pixels used
at each gray level, by placing the vertical lines in the right hand
distribution plots at equally spaced quantiles of the matrix values
distribution, as shown in the right hand panel of Figure \ref{figF.B3:ToyHeatMapScalingQuantile}.
Note that the corresponding heat-map on the left now strongly reveals
the regular background pattern. As with many choices among OODA methods,
there is a trade-off of which one should be aware: the clear view
of the background comes at the cost of much less contrast between
the heights of the major peak and the three other large ones.

\begin{figure}[H]
\includegraphics[scale=0.9]{ChapterF/OODAbookChpFFigB3ToyHeatMapQuantile}

\caption{Another heat-map visualization of the same data matrix from Figures
\ref{figF.B1:ToyHeatMapScalingLinear} and \ref{figF.B2:ToyHeatMapScalingLog}.
Highlights different structure in the data using the even spread of
the gray scale colors provided by quantile scaling. \label{figF.B3:ToyHeatMapScalingQuantile}}
\end{figure}
Notice that quantile scaling is the same whether applied to the original
data or the log transformed version, because the log is a monotone
transformation. However the matrix value distributions in the right
panels above are much easier to understand on the log scale (used
in Figures \ref{figF.B2:ToyHeatMapScalingLog} and \ref{figF.B3:ToyHeatMapScalingQuantile})
than on the linear scale of Figure \ref{figF.B1:ToyHeatMapScalingLinear}.

When choosing a heat-map visualization tool, at least having the option
for some type of visualization of these heat-map gray level (or color)
distributions can be very useful. An attractive option is a histogram,
where each bar represents a gray level (or color) with each color
used for its corresponding bar. Best insight comes from taking the
area of each bar to be proportional to the number of pixels using
that color. This is the same data display principle that underlies
density estimation\index{density estimation} as discussed in Section
\ref{secN-A:NoHisto}. Such visualizations of color distribution appear
in the lower right panels of Figures \ref{figF.F1:ComboViewInputData}
- \ref{figF.F4:ComboViewMode2}.

For data matrices taking both positive and negative values where it
is important to highlight 0, a common two color choice codes 0 as
black, with shades of red and green indicating magnitude and direction
of departures from 0. This color scheme has the disadvantage of being
inaccessible to red-green \index{color blind}color blind people,
which is nearly 10\% of the male American and European populations.
An alternate choice illustrated in Figure \ref{figF.C1:ToyHeatMap},
codes 0 as white with shades of red and blue for magnitudes and directions.
In some fields (e.g. economics or finance) red is frequently used
for negative values (suggesting financial loss), while in others (e.g.
climatology) it is more natural to use red for positive (corresponding
to hotter temperatures). In Figure \ref{figF.C1:ToyHeatMap} the former
scheme is used. Note that the actual coded values are shown using
a color bar on the right side of the heat-map.

The main point of Figure \ref{figF.C1:ToyHeatMap} is to show that
while a heat map can provide very useful insights in some situations,
it can miss important structure in others. This is done using the
Two Class Gaussian data\index{Two Class Gaussian data|(} set, which
is a matrix with $n=200$ data object vectors of length $d=20,000$,
which is common in bioinformatics. To enhance visible structure as
demonstrated in Figure \ref{figF.A:ToyHeatMapCluster}, hierarchical
clustering\index{hierarchical clustering|)}\index{clustering} using
average linkage\index{linkage!average} and Euclidean distance (a
common choice for visualization in bioinformatics as noted above)
has been applied to both rows and columns. A challenge to heat map
visualization of data sets of this size is that the pixel capability
of most types of displays tend to be in the very low thousands, which
is very inadequate in this case. Hence, most heat-map visualization
software packages provide convenient facilities for zooming in on
appropriate parts of the display. Here that is done by simply showing
the first 200 rows of this much longer matrix. Other subsets show
similar apparent pure noise.

\begin{figure}[H]
\includegraphics[scale=0.9]{ChapterF/OODAbookChpFFigC1HeatMap}

\caption{Heat-map image view of Two Class Gaussian data set, taking on both
positive and negative values. White indicates matrix entries that
are 0, with blue and red for positive and negative values, respectively.
It is hard to see much structure, even though rows and columns have
been clustered. \label{figF.C1:ToyHeatMap}}
\end{figure}
As suggested by the name, there actually is structure in this Two
Class Gaussian data matrix. In particular it was generated as the
first $n_{1}=100$ columns having all entries distributed as independent
$N(0.04,1)$, with the next $n_{2}=100$ columns drawn independently
from the $N(-0.04,1)$ distribution. While this perhaps small true
underlying signal is dominated by the noise in Figure \ref{figF.C1:ToyHeatMap},
it can be seen in other types of visualizations, such as the PCA scatterplot
matrix\index{scatterplot matrix|(} as shown in Figure \ref{figF.C2:ToyPCA}.
Here the two subpopulations (labeled as magenta circles and green
plus signs) are visually distinct, showing that the PCA scatterplot
matrix view can reveal structure that is not apparent in the heat-map
image. 

\begin{figure}[H]
\includegraphics[scale=0.85]{ChapterF/OODAbookChpFFigC2PCA}

\caption{PCA scatterplot matrix view of Two Class Gaussian data\index{Two Class Gaussian data|)}
set in Figure \ref{figF.C1:ToyHeatMap}. Two subpopulations, indicated
with different colors and symbols, are visibly quite distinct, despite
the high noise level, in sharp contrast with the heat-map view. \label{figF.C2:ToyPCA} }
\end{figure}
However, as seen in Section \ref{secB3:Confirmatory-Analysis}, one
should be skeptical as to whether such visually apparent aspects of
a scatterplot represent true underlying structure, or are merely irreproducible
sampling artifacts. In this case, the DiProPerm\index{DiProPerm}
analysis of this Two Class Gaussian data set in Figure \ref{figL.A:ToyPCADiProPerm}
shows the difference between these groups is strongly significant
(p-value much less than 0.001), as it should because there is a systematic
difference between the subgroups.

In summary the heat-map visualization tool has proven to be very powerful
in the hands of experienced analysts. As discussed at the beginning
of Section \ref{secF.A:Heatmap-Views}, it has a strong track record
of finding clusters and their drivers in a single visualization. But
as seen in Figure \ref{figF.C2:ToyPCA}, it can also miss important
structure in high noise situations. One more point to keep in mind
is that apparent clusters in a heat-map may not represent reproducible
discoveries. A hypothesis test that is well suited for analyzing the
statistical significance of clusters (especially in high dimensional
contexts) is SigClust\index{SigClust}, discussed in Section \ref{secL2:SigClust}.

\section{Curve Views of Matrices and Modes of Variation\label{secF.B:Curve-Views}}

Curve views\index{data matrix!curve view} (standard in FDA\index{functional data analysis})
reveal quite different aspects of variation in a sample of data objects,
as seen in many situations in earlier chapters. The value of such
views includes a clear display of modes of variation\index{modes of variation|(},
as defined in Section \ref{subsecBA.1.1:FormalizingModes}. Some specific
examples follow.
\begin{itemize}
\item \emph{Figure \ref{fig:Mortality-Raw}, Spanish Mortality data}.\index{Spanish Mortality data}
These displays demonstrate the value of the log transformation\index{transformation},
in terms of making much of the variation in the data available to
the analysis.
\item \emph{Figure \ref{fig:Mortality-Colored}, Spanish Mortality data}.
This curve view graphic, together with the previous figure highlights
the value of careful color coding for clear insights about (for example)
time structure in the data curves. 
\item \emph{Figure \ref{fig:Mortality-Mean-Resid}, Spanish Mortality data}.
Here the focus is on variation about the mean, revealing that many
of the age impacts on mortality have been constant over time. 
\item \emph{Figure \ref{fig:Mortality-PC1}, Spanish Mortality data}. The
main mode of variation is an overall improvement in mortality, which
is most dramatic for the young, together with an improvement in the
keeping of records.
\item \emph{Figure \ref{fig:Mortality-PC2}, Spanish Mortality data}. The
second mode of variation is a contrast between the 20-45 year old
males with the rest of the population. Recall this mode reflected
pandemic, war and automobile effects.
\item \emph{Figure \ref{fig:TIC-spectra}, TIC spectra}\index{TIC Curves data}.
This curve view display shows that alignment can be a challenging
problem, that is very effectively addressed by the Fisher-Rao approach
to curve registration\index{curve registration!Fisher-Rao}.
\item \emph{Figure \ref{fig:Sim-Warps}, Bimodal Phase Shift data}\index{Bimodal Phase Shift data}.
These curves demonstrate how appropriately chosen warping\index{warping}
functions can be used to align peaks in curve registration, and to
visually reflect phase variation (lower right panel).
\item \emph{Figure \ref{fig:10d-ToyPCA}},\emph{ Tilted Parabolas data}.\index{Tilted Parabolas data}
There are a number of curve displays in this analysis which give a
clear indication of the variation in the full data, the mean residuals,
and several principal components together with their respective residuals.
The left plot in the second row shows that vertical shift is the dominant
mode of variation (a common aspect of many functional data sets that
is more deeply explored in Section \ref{subsecP-A-A:Data-Centering}).
That on the third row reveals the less obvious tilting mode of variation.
The bottom row suggests that the remaining variation is random.
\item \emph{Figure \ref{fig:50d-ToyPCA}, Twin Arches data}\index{Twin Arches data|(}.
Here the curve view shows two non-obvious modes of variation, with
the first component reflecting a peak-valley versus valley-peak mode,
while the second mode is peak-peak versus valley-valley.
\item \emph{Figures \ref{fig:LungCancerLogCurves} and }\ref{fig:LungCancerLogCurvesBrushed}\emph{,
\index{Lung Cancer data}Lung Cancer data}. The first shows a large
amount of variation that is hard to visually parse. The cluster based
\index{brushing}brushing (i.e. coloring) of these curves shown in
the second reveals important modes of variation as discussed in Section
\ref{secB2:Data-Visualization}. 
\item \emph{Figures \ref{fig:ToyPCAnScaling} and \ref{fig: ToyPCAnScalingCorr},
Two Scale Curves data\index{Two Scale Curves data}}. These views
of the data demonstrate the potential impact of data normalization\index{normalization}
(through standard deviation scaling, as done in unit free correlation
matrix approaches) on PCA. In particular this choice can highlight
completely different modes of variation.
\item \emph{Figure \ref{figD-M:ToyShiftsnFRReg}}. Here studying curves
that are mostly horizontally shifted clarifies the impact of peak
alignment on the pointwise mean of a bundle of curve data objects.
\end{itemize}
These examples demonstrate the potential benefits of curve views when
the data objects are curves, especially when used in tandem with PCA,
or some other approach to revealing modes of variation. However as
with any OODA method, there are situations where this approach is
less useful, for example when there is a large amount of noise present.
This is illustrated using the Two Clusters toy data (shown in the
right panel of Figure \ref{figF.A:ToyHeatMapCluster} where both rows
and columns have been hierarchically clustered) by a curve view in
Figure \ref{figF.D:ToyPoorCurveView}. Following the convention in
this book, columns of the matrix are the (curve) data objects, so
the horizontal axes in Figure \ref{figF.D:ToyPoorCurveView} show
the indices of the rows. The structure visible in the heat-map in
the right panel of Figure \ref{figF.A:ToyHeatMapCluster} can also
be seen in Figure \ref{figF.D:ToyPoorCurveView}, but the interpretation
takes more effort. In particular, the lighter top half of the former
appears as higher values of the mean curve on the left shown in the
top center panel of the curve view. The brighter gray levels on the
right half of the heat-map (i.e. one cluster of the column data objects)
appear as the upper set of cyan curves in the lower left panel of
Figure \ref{figF.D:ToyPoorCurveView}, and as the cyan cluster in
the scores shown in the lower right panel. The lower central panel
shows the mean plus and minus the extreme curves, which contain all
of the structure available in the heat-map. The top right panel shows
the scree plot (introduced in Figure \ref{figB-Y:Screeplots}), revealing
that the noise level is indeed quite high in this example. While both
views distinguish the underlying structure in the data set, in this
case the heat-map is much easier to interpret.

\begin{figure}[H]
\includegraphics[scale=0.9]{ChapterF/OODAbookChpFFigDToyPoorCurveViewCMYK}

\caption{Curve view of the Two Clusters data set from Figure \ref{figF.A:ToyHeatMapCluster}.
Shows that curve view may not be as useful as a heat-map view in high
noise contexts. \label{figF.D:ToyPoorCurveView}}
\end{figure}
In very high noise data sets, such as the Two Clusters data\index{Two Clusters data}
in Figures \ref{figF.A:ToyHeatMapCluster} and \ref{figF.D:ToyPoorCurveView},
the noise level can have a substantial impact on PCA, in terms of
both loadings and scores. This motivated a large amount of research
in the FDA\index{functional data analysis} community on integrating
smoothing\index{smoothing} with analyses such as PCA. See \citet{RamsaySilverman2002,RamsaySilverman2005}
and \citet{FerratyVieu2006} for good overviews. Additional aspects
are treated in \citet{ZhangJT2013ANOVA4FDA} and \citet{Kokoszka2013FDAbook}.
A particularly strong theoretical treatment of FDA can be found in
\citet{HsingEubank2015TheoryFDA}. A simple approach is to simply
smooth the data object curves before doing PCA, see Chapter \ref{chapN:Smoothing}
for further discussion of smoothing methods. One could also do PCA,
and then smooth the loadings curves. However generally better is to
integrate the smoothing with the PCA, which can be done by various
methods as detailed in the above references. Integration of smoothing
methods with PCA is also very useful for functional data sets with
curves evaluated at sparsely sampled or irregular points. An important
method in such situations is the PACE algorithm introduced in \citet{YaoMuellerWang2005PACE}.
\citet{DelaigleHall2016FragmentedFDA} propose an approach to handling
a surprisingly high level of sparse sampling.

Less well understood is that in standard low to moderate noise situations,
there is not much gained from explicitly integrating smoothing methods
with PCA, because simple PCA frequently has an implicit smoothing
effect. This can be seen in the Tilted Parabolas data\index{Tilted Parabolas data}
in Figures \ref{fig:10d-ToyPCA} and \ref{fig:10d-ToyPCA-MinWhiteSpace},
and in the Twin Arches data in Figure \ref{fig:50d-ToyPCA}, where
the substantial level of noise visible in the raw data appears in
the residual, not in the main components. This appears to be caused
by the noise being mostly averaged out when the PCA calculation estimates
the smooth underlying loadings vectors.

\section{Data Centering and Combined Views\label{secF.E:ComboView-Centering}}

The \emph{mean centering}\index{mean!center|(} of data is an issue
that is very frequently dismissed as trivial. The impact of centering\index{centering|(}
on an OODA seems fairly obvious and routine from a curve viewpoint,
as illustrated for example in Figure \ref{fig:Mortality-Mean-Resid}.
However, when contemplating heat-map views the effects can be surprisingly
important and even challenging to intuitively comprehend. In particular,
that view motivates the consideration of both horizontal and vertical
mean centering of the data matrix. The intuitive impact of these can
be surprisingly elusive. That issue has been more deeply investigated
in \citet{Prothero2021Centering}.

The importance of doing a centering operation (before looking for
modes of variation), from the curves as data objects viewpoint (e.g.
in FDA), is illustrated in Figure \ref{figP-B:SVDvsPCA}. This contrasts
PCA with a fully uncentered analysis based on a direct Singular Value
Decomposition\index{Singular Value Decomposition} (SVD) of the data
matrix. Many aspects of SVD and its relation to PCA are studied in
detail in Section \ref{subsecP-A-B:SVD}. 

Figure \ref{figP-B:SVDvsPCA} provides a comparison between uncentered
SVD and PCA in the context of the simple 2-d Toy data example from
Section \ref{subsecBA.A.B:ObjFeatExample}. Using the terminology
of Section \ref{secBA.1:Foundations_Terminology}, the data are shown
in a feature space\index{feature space} view as the black circles
in the top two panels. The top left panel is a view of the SVD approximation
of this data set. The red line is the one dimensional subspace (i.e.
line through the origin) that best fits the data. The magenta plus
signs are the projection\index{projection|(}s of the data onto that
subspace, and are the best rank one approximation of the data. In
particular this direction minimizes the sum of squares of the projected
residuals (shown as cyan lines). The first SVD \index{scores|(}scores
are the coefficients of these projections, which appear on the horizontal
axis in the left lower panel. The signed lengths of the cyan lines
are also the coefficients of the projections onto the second singular
vector, i.e. the second SVD scores, which are used on the vertical
axis in the lower left panel. Because SVD ignores the center of the
data, it fails to efficiently summarize and display the dominant mode
of variation in this data set (i.e. major axis of the ellipse of data
which is essentially vertical). Instead that notion of main variation
is split between both of these modes.

The right hand panels of Figure \ref{figP-B:SVDvsPCA} illustrate
the PCA of the same 2-d Toy data\index{2-d Toy data} set. Recall
the difference with SVD is the centering of the analysis at the mean
data object (shown as the green plus). This results in the best approximating
line, shown in red, being now chosen from direction vectors based
at that point. The rank 1 PCA approximations are shown as magenta
plus signs, which clearly provide a much better summarization of the
data than is provided by SVD, because this direction maximizes the
variance of the projections. This direction now appropriately reflects
the dominant vertical mode of variation. In particular, these rank
one approximations now lie in the middle of the point cloud. This
demonstrates a clear value to data object mean centering. The coefficients
of these projections are the PC1 scores, plotted on the horizontal
axis of the scores scatterplot in the bottom right panel. The cyan
lines in the top right panel show the residuals of this approximation
(which of course have minimal sum of squares). The signed lengths
of these are the PC2 scores used on the vertical axis in the bottom
right panel. 

\begin{figure}[H]
\includegraphics[scale=0.85]{ChapterP/OODAbookChpPFigBSVDvsPCA}

\caption{Toy example comparing SVD and PCA for the 2-d Toy data. Top panels
are feature space views of the data objects (black circles) together
with the SVD (left) and PCA (right) approximations. Corresponding
scores representations appear in the bottom panels. Shows the column
mean object centering implicit to PCA allows better low rank approximation
of the data, and also uncorrelated scores scatterplots.\label{figP-B:SVDvsPCA}}
\end{figure}
Figure \ref{figP-B:SVDvsPCA} also highlights another aspect of PCA
(relative to SVD of uncentered data) that is worth noting: the scores
are uncorrelated. Essentially this is because when there is correlation,
as shown in the lower left panel, the variance of the larger projection
can be made even larger by appropriate rotation. Similarly the sum
of squared residuals (shown as cyan in the top row) can be made smaller
by this rotation. It is perhaps surprising, and generally not widely
acknowledged, that the lack of correlation in PC scores is a consequence
of the column data object mean centering. The mathematics behind these
observations is given in the calculation (\ref{eqP.AE:OrthoScores})
in Section \ref{subsecP-A-B:SVD}.

There are several reasonable notions of centering (i.e. subtracting
some type of mean) of a data matrix $\boldsymbol{X}$ (using notation
from (\ref{eqB.E:Matrices})). From the default OODA perspective of
column vectors as \index{data objects!column vectors}data objects,
a natural definition of center is the \emph{mean vector} 
\begin{equation}
\overline{\boldsymbol{x}}_{CO}=\left(\begin{array}{c}
\overline{x}_{1,A}\\
\vdots\\
\overline{x}_{d,A}
\end{array}\right)\in\mathbb{R}^{d},\label{eqF1:ColumnObjectMean}
\end{equation}
where $\overline{x}_{i,A}=n^{-1}\sum_{j=1}^{n}x_{i,j}$ for $i=1,...,d$.
A potential ambiguity (which tends to generate confusion in discussions
about data analysis) is whether $\overline{\boldsymbol{x}}_{CO}$
should be called the ``column mean'', as it is the mean vector of
the column vectors, or should be called the ``row mean'' because
it is the vector whose entries are the means of the entries in the
corresponding rows of the data matrix $\boldsymbol{X}$. The OODA
resolution of this problem is to use the name \emph{column object
mean\index{column object mean|(}} for $\overline{\boldsymbol{x}}_{CO}$,
hence the subscript of $CO$. This choice of terminology is deliberately
intended to keep the focus of the discussion on the data object column
vectors, and not on individual matrix entries. 

While column data object vectors are the main focus of OODA, the row
vectors that make up the data matrix $\boldsymbol{X}$ are sometimes
also of interest, e.g. when considering multi-block methods as in
Section \ref{secP-B:TwoBlockDecomp}. Using terminology from Table
\ref{tab2-A:Synonyms} OODA terminology for the rows of $\boldsymbol{X}$
is \emph{\index{row trait!vector|(}row trait vectors}. Here the more
biological synonym of ``trait\index{trait|(}'' instead of ``feature\index{feature}''
is used because the term ``feature vector'' already has the common
meaning of ``data object column vector''. The term ``trait vectors''
for row vectors of the data matrix avoids that ambiguity. Averaging
these row vectors gives the \emph{row trait mean\index{row trait!mean|(}
vector}, 
\begin{equation}
\overline{\boldsymbol{x}}_{RT}^{t}=\left(\begin{array}{ccc}
\overline{x}_{A,1} & \cdots & \overline{x}_{A,n}\end{array}\right),\;\overline{\boldsymbol{x}}_{RT}\in\mathbb{R}^{n},\label{eqF2:RowFeatureMean}
\end{equation}
which has entries $\bar{x}_{A,j}=d^{-1}\sum_{i=1}^{d}x_{i,j}$ for
$j=1,\cdots,n$. 

As noted above deep study of these notions of center, including their
representations as projections, can be found in Section \ref{subsecP-A-A:Data-Centering}.
The rest of this section considers how centering by subtracting one
or the other (or both) of the column object mean (\ref{eqF1:ColumnObjectMean})
and the row trait mean (\ref{eqF2:RowFeatureMean}) impacts a data
matrix.

Relationships between heat-map and curve views of a common data set
are explored using a \emph{combined view\index{combined view|(}\index{data matrix!combined view}}
in Figure \ref{figF.F1:ComboViewInputData}. This is an analysis of
the $d=50$ dimensional Twin Arches data set that was studied in Figures
\ref{fig:50d-ToyPCA} and \ref{fig:50d-ToyPCAScatMat}. That data
set was generated with theoretical mean 0 in each coordinate, but
due to sampling variation the column object sample mean vector $\overline{\boldsymbol{x}}_{CO}$,
shown in the top center panel of Figure \ref{fig:50d-ToyPCA}, is
only nearly but not exactly equal to $\boldsymbol{0}_{d,1}$. Various
views of the input data set are shown in Figure \ref{figF.F1:ComboViewInputData}.
The upper right panel is essentially the same as the upper left panel
of Figure \ref{fig:50d-ToyPCA} except it has been transposed (so
the coordinate axis is now vertical, with curve height on the horizontal
axis), and the curves are differently colored. These same curves (actually
just vectors) are grouped into a matrix (with corresponding vertical
axes) in the heat-map view in the top left panel, using intensities
of red (blue) to indicate negative (positive respectively) curve heights
with white for 0. The four clusters revealed in the middle top panel
of Figure \ref{fig:50d-ToyPCAScatMat} appear as four distinct vertical
red-blue color patterns in the columns of the heat-map. Recall from
the middle row of Figure \ref{fig:50d-ToyPCA} that two large clusters
had arches that followed either a valley-peak or a peak-valley pattern.
These appear in the right hand part of the heat-map, which make it
clear that these clusters drove the first PC because there are more
column data objects in them. Similarly the last row of Figure \ref{fig:50d-ToyPCA}
indicates curves with peak-peak and valley-valley patterns, which
appear in the first two smaller clusters of column objects on the
left of the heat map. The color scheme of the curves in the top right
panel follows the rainbow color bar in Figure \ref{fig:Mortality-Colored},
from magenta through blue, green and yellow to red, in the order of
the columns of the matrix. The yellow-red curves are plotted last
and follow the valley-peak pattern as shown in the far right cluster
of column objects in the heat-map. Second to last are the cyan-green
curves with a peak-valley pattern that is similarly consistent with
the third set of columns in the top left heat map. The other colors
only appear intermittently, as those curves are mostly overplotted
by the latter two clusters. 

\begin{figure}[H]
\includegraphics[scale=0.85]{ChapterF/OODAbookChpFFigF1ComboViewInputDataCMYK}

\caption{Combined view of the Twin Arches toy data set from Figures \ref{fig:50d-ToyPCA}
and \ref{fig:50d-ToyPCAScatMat}. Upper left is a heat-map image view
(red for negative, blue for positive values). Upper right is the columns
as data object curves view (using a rainbow color scheme), transposed
for good correspondence with the heat map. Lower left shows the corresponding
view of the row trait vectors as curves (with a heat ordered color
palette). Lower right is distribution of colors in the heat-map. \label{figF.F1:ComboViewInputData}}
\end{figure}
The bottom left panel of Figure \ref{figF.F1:ComboViewInputData}
shows the same data using a parallel coordinates plot of the row trait
vectors. Again the horizontal axis is deliberately chosen to correspond
to the horizontal axis of the heat map. To make clear the contrast
with other panels, a completely different color scheme is employed,
which is heat ordered from black through red and orange to light yellow.
A key point is that the common variation in the heat map is expressed
in quite different ways, although there is an important duality between
them. In particular, instead of having four clusters determined by
just two peak-valley patterns (as for the column data objects discussed
above), from this viewpoint the upper left heat-map shows that the
row trait vectors naturally fall into only two clusters determined
by more complicated peak-valley patterns. The cluster of rows appearing
in the bottom of the heat-map follows a peak-valley-peak blue-red
pattern whose curves are highlighted in the bottom left using cooler
black and dark red colors. The peak-valley-peak pattern in the curve
heights reflects that which is shown with blue-red colors in the upper
left panel. The other cluster follows a different peak-valley-peak-valley
pattern (again both in the top row cluster in the heat-map as well
as in the curves) shown using curves colored with hotter orange and
yellow shades in the bottom left panel. The duality between the column
object and row trait curve representations is that clusters in one
correspond to curve structure in the other. 

The bottom right panel of Figure \ref{figF.F1:ComboViewInputData}
shows the distribution of red-white-blue colors used for the pixels
in the heatmap. For each color, the height of the bar (and hence its
area) indicates the frequency of pixels in the heatmap of that color.
The use of this type of graphic was discussed in Section \ref{secF.A:Heatmap-Views}.

The toy example of how SVD compares to PCA in Figure \ref{figP-B:SVDvsPCA}
demonstrated the benefits of centering the data by subtracting the
column object mean $\overline{\boldsymbol{x}}_{CO}$ (i.e. the mean
of the column vectors). However the heat map view of the data begs
the question of why not also center the other way by subtracting the
mean of the row trait vectors $\overline{\boldsymbol{x}}_{RT}^{t}$?
Note this is the same as subtracting the mean of the entries from
each column. Data adjustments of this type are particularly important
for various types of bioinformatics data. For example most gene expression
measurements involve amplification of quantities such as RNA which
essentially involves an uncontrollable random multiplicative factor
for each case. When such data are analyzed on a log scale\index{scale}
(the usefulness of this is demonstrated in Figures \ref{fig:LungCancerLogCurves}-\ref{fig:LungCancerLogCurvesBrushed}),
such effects can be mitigated by row trait mean centering, which essentially
controls for that random amount of scaling in the original measurements. 

The residuals from doing both types of centering, called \emph{double
\index{centering!double|(}centering} in Section \ref{subsecP-A-A:Data-Centering},
for the toy data in Figure \ref{figF.F1:ComboViewInputData} are shown
in Figure \ref{figF.F2:ComboViewOAMeanResid}. In this case, the column
object centering had a negligible impact because as noted above $\overline{\boldsymbol{x}}_{CO}\approx\boldsymbol{0}$,
however the row trait centering does change the visualization. Note
that the column object curves view of these mean residuals in the
upper right panel is quite different from the column object mean residuals
shown in the upper right panel of Figure \ref{figF.F1:ComboViewInputData}.
The predominantly cyan-green and orange-red bundles of curves are
the same as in the upper left of Figure \ref{figF.F1:ComboViewInputData},
because the means of the entries of those column vectors are all essentially
0 (apparent in both top panels of Figure \ref{figF.F1:ComboViewInputData}).
The two new bundles that were previously hidden under those bundles
(due to being plotted earlier in the process of constructing the upper
right panel of Figure \ref{figF.F1:ComboViewInputData}), but have
been pulled out by the row trait centering since the mean of their
entries is not 0. In particular, each entry of the row trait mean
vector $\overline{\boldsymbol{x}}_{RT}^{t}$ from (\ref{eqF2:RowFeatureMean})
is calculated as the average of the entries of the corresponding heat
map column. For the first (left-most) cluster of columns in the heat
map, the matrix entries are all between 0 and 10 (apparent in both
the upper left heat map and the lower left row curve view of Figure
\ref{figF.F1:ComboViewInputData}) so the averages are around 5. The
impact of subtracting 5 from each column of the matrix moves the color
range {[}0,10{]} to {[}-5,5{]} in the upper left heat map, and similarly
moves the bundle of curves in the lower left panel downwards, while
in the upper right panel, this shift of 5 moves only the predominantly
magenta bundle of curves to the left, i.e. out from under the other
bundles. Similarly, the second cluster in the heat map has column
averages around -5, so that cluster has opposite shifts in each panel.
In particular in the heat map, the color range shifts from {[}-10,0{]}
to {[}-5,5{]}, with a corresponding shift in the lower left panel,
as well as pulling the predominantly blue second cluster out from
the other bundles towards the right in the upper right panel. A key
lesson is that while column object mean centering results in a rigid
shift of the full set of column object curves, the row trait mean
centering operates on individual column object curves.

\begin{figure}[H]
\includegraphics[scale=0.85]{ChapterF/OODAbookChpFFigF2ComboViewOAMeanResidCMYK} 

\caption{Combined view of Twin Arches data, after both column object and row
trait centering (i.e. double centering). Shows interpretation of both
types of centering is not always intuitive.\label{figF.F2:ComboViewOAMeanResid} }
\end{figure}
As will be seen in Section \ref{subsecP-A-B:SVD}, for double centered
data (both row and column centered), computation of the PCA modes
of variation for \index{row trait!vector|)}row trait data objects
and for column object data objects are very closely linked. Essentially
the loadings\index{loadings} for one provide the scores for the other.
The first of these PCA modes, for the Twin Arches data of Figures
\ref{figF.F1:ComboViewInputData} and \ref{figF.F2:ComboViewOAMeanResid}
is studied using a similar combined view in Figure \ref{figF.F3:ComboViewMode1}.
The column object view in the upper right is essentially the transpose
of the middle left first mode of variation shown in Figure \ref{fig:50d-ToyPCA},
although this color scheme is more interpretable. The reason for this
is clear from either panel on the left side. In particular this mode
is strongly driven by points in the third and fourth column clusters
of the heat map image and essentially ignores the first two clusters
which were strongly affected by the row trait mean adjustment. The
top left heat map image also illustrates how PCA can be viewed as
the best rank one approximation. In particular, the shown matrix explains
as much variation in the full set as possible with a rank one matrix.
Note that the larger number of nearly white pixels (on the left side
of the upper right heat-map) is clearly displayed as nearly zero values
in the color distribution plot in the lower right. Also note the much
higher vertical scale compared to the color distribution (lower right
panel) in Figure \ref{figF.F2:ComboViewOAMeanResid}.

\begin{figure}[H]
\includegraphics[scale=0.85]{ChapterF/OODAbookChpFFigF3ComboViewMode1CMYK}

\caption{First mode of variation (best rank one approximation) based on double
mean centered Twin Arches data from Figure \ref{figF.F2:ComboViewOAMeanResid}.
Upper right panel shows the column object mode of variation is similar
to that in Figure \ref{fig:50d-ToyPCA}. Left two panels reveal this
is driven by the large two column object clusters revealed above.
\label{figF.F3:ComboViewMode1}}
\end{figure}
The combined view of the second mode of variation based on double
centering is shown in Figure \ref{figF.F4:ComboViewMode2}. In contrast
to the first mode shown in Figure \ref{figF.F3:ComboViewMode1}, the
column object view in the upper right panel appears to be quite different
from the mode displayed in the lower left panel of Figure \ref{fig:50d-ToyPCA}.
But both panels on the left make it clear why. They both reveal that
this variation appears in nearly only the first two column object
clusters, i.e. in the curves numbered 1-16 on the horizontal axes.
At first glance the curves in the top right panel appear to exhibit
twice the frequency of those in the lower left of Figure \ref{fig:50d-ToyPCA},
but a careful inspection reveals that this is just the result of separately
row trait mean centering each of those two clusters. From the matrix
approximation viewpoint, this rank one matrix is the product of a
row and column vector, each of which is orthogonal to the corresponding
vector in Figure \ref{figF.F3:ComboViewMode1}, which best approximates
the remaining variation. Note that both axes of the color distribution
in the lower right panel are quite different from the preceding cases
showing even more nearly white pixel\index{pixel|)}s.

\begin{figure}[H]
\includegraphics[scale=0.85]{ChapterF/OODAbookChpFFigF4ComboViewMode2CMYK}

\caption{Second mode of variation (rank one approximation) based on double
mean centered Twin Arches data from Figure \ref{figF.F2:ComboViewOAMeanResid}.
The column object mode of variation shown in the upper right panel
this time appears different from that in Figure \ref{fig:50d-ToyPCA}.
Left two panels again explain this via insights about clustering of
column objects. \label{figF.F4:ComboViewMode2}}
\end{figure}
In\index{combined view} the examples shown in Figures \ref{figF.F1:ComboViewInputData}
- \ref{figF.F4:ComboViewMode2}, basing the analysis on double mean
centering gives a very insightful decomposition of the structure of
the data set. However, this is not always the case, as illustrated
in Figure \ref{figG.G:LungCancerCentered}. That studies the same
Lung Cancer RNAseq\index{RNAseq data} gene expression curves shown
in Figures \ref{fig:LungCancerLogCurves} - \ref{fig:LungCancerLogCurvesBrushed}.
The left panel of Figure \ref{figG.G:LungCancerCentered} shows the
column object mean\index{column object mean|)} centered version of
the gene expression curves, using the colors from Figure \ref{fig:LungCancerLogCurvesBrushed}.
Note that the essentially unexpressed cases colored red, are still
clearly distinguished. Furthermore the critical gold versus blue distinction
in the exon near the right side is also still very clear. However,
there is some loss of insight relative to the very clear explanation
available from the Figure \ref{fig:LungCancerLogCurvesBrushed} view.
This loss of interpretability becomes much worse in the double centered
view shown in the right panel of Figure \ref{figG.G:LungCancerCentered}.
This time the impact of row trait mean centering (with strong potential
for moving clusters as demonstrated in the top right panel of Figure
\ref{figF.F4:ComboViewMode2}) moves these clusters in a way that
strongly impairs seeing the key lessons in this data set. In particular
the fact that the (unexpressed) red cases are substantially lower
than the others mostly disappears because each curve has separately
had the mean of its entries subtracted which obscures important aspects
of the data.

\begin{figure}[H]
\includegraphics[scale=0.7]{ChapterF/OODAbookChpFFigGlungCancerCentered}

\caption{Centered views of Lung Cancer data from Figure \ref{fig:LungCancerLogCurvesBrushed}.
Column object centering on the left and double centering\index{centering!double|)}
on the right. Shows how centering can give less interpretable views
of data. \label{figG.G:LungCancerCentered}}
\end{figure}
An important lesson from these examples is that substantial thought
should be devoted to centering of data. Both column object centering
and row trait centering have the potential to either reveal or obscure
important insights. Further centering\index{centering|)} issues,
which include various projection viewpoints, are studied in greater
depth in Section \ref{subsecP-A-A:Data-Centering}.\index{mean!center|)} 

\section{Scatterplot Matrix Views of Scores\label{secF.C:Scatterplot-Matrix-Views}}

Scatterplot matrix views provide yet another frequently useful type
of insight, this time highlighting the \emph{relationships between\index{data objects!relationships between}
data objects}. This has also been seen many times in the above examples
as detailed here.
\begin{itemize}
\item \emph{Figure \ref{fig:50d-ToyPCAScatMat}, Twin Arches example.} This
view provided clear presentation of the important structure in the
data set that was challenging to understand from the curve views of
this data in Figure \ref{fig:50d-ToyPCA}. In particular, it revealed
that there are four clusters in the data, whose 1-d projections generated
the perhaps surprising effects shown there. That example also illustrated
axis scaling issues, where the PC3 and 4 axes were rescaled to minimize
white space, at the cost of not clearly revealing that those scores
were an order of magnitude smaller than the others, i.e. the noise
component was generated from an isotropic Gaussian distribution.
\item \emph{Figure \ref{fig:LungCancerPCA}, Lung Cancer data}. This PCA
scatterplot matrix revealed the important structure of clusters in
the data, which motivated the development of a novel method for finding
alternate splicing using gene expression data.
\item \emph{Figure \ref{fig:LungCancerPCABrushed}, Lung Cancer data}\index{Lung Cancer data}.
This brushed\index{brushing} version of Figure \ref{fig:LungCancerPCA}
defined a color scheme whose application led to Figure \ref{fig:LungCancerLogCurvesBrushed},
which revealed the nature of the clusters in this data set, in terms
of alternate splicing. This led to the invention of the SigFuge method,
that was effective at finding new alternate splicing events.
\item \emph{Figure \ref{figB-Q:PanCanPCA}, Pan Cancer data}\index{Pan Cancer data}.
Here is a PCA scatterplot aimed at showing the limitations of PCA,
using a cancer gene expression data set. Colors and symbols indicate
cancer types, and reveal some differences, although also substantial
overlap of some of the classes. 
\item \emph{Figure \ref{figB-R:PanCanDWD}, Pan Cancer data}. This demonstrated
how computing scores by projecting the data onto non-PCA directions
can give much better visual contrast between cancer types. This used
DWD\index{Distance Weighted Discrimination|(} directions (discussed
in Chapter \ref{chapJ:Classification}) trained on pairs of individual
cancer classes, which directly target maximal separation of the classes.
\item \emph{Figure \ref{fig:HiDimToyPCA }, High dimensional Gaussian example}.
Here is another example designed to reveal limitations of PCA. No
apparent structure is apparent in this PCA scores scatterplot, especially
in terms of comparison of the two classes.
\item \emph{Figure \ref{fig:HiDimToyDWD}, High dimensional Gaussian example}.
This again shows that the DWD direction, which deliberately optimizes
for separation between classes, reveals a strong difference between
them. Such views generally improve on PCA for understanding potential
class differences. Note that the DWD scores exhibit substantially
less variation than is apparent in the PC directions. This is again
an axis scaling issue, that could be highlighted at the cost of introducing
white space.
\item \emph{Figure \ref{figD-A:DrugPCARaw}, Drug Discovery data\index{Drug Discovery data}}.
The PCA view revealed that the raw data were dominated by outliers,
which distracted from the important activity differences which were
the main motivation of the analysis.
\item \emph{Figure \ref{figD-J:DrugVDBinaryPCA}, Drug Discovery data}.
Here the impact on the PCA view of data cleaning and focusing only
on the binary variables is revealed. Shows much richer structure in
the data, and reveals interesting insights about activity differences,
shown with color and symbols.
\item \emph{Figure \ref{figD-K:DrugVDnonBinaryPCA}, Drug Discovery data}.
This shows that the appropriately cleaned and transformed non-binary
variables do a much better job of understanding activity differences.
\item \emph{Figure \ref{figF.C2:ToyPCA}, Two Cluster example}. Here it
is seen how scatterplot matrices can reveal structure in data that
is not easily available from heat-map\index{heat-map|)} views.
\end{itemize}
In most of the above scatterplot matrices, the graphics are straightforward
because the directions that determine the axes are orthogonal. But
this is not always the case, for example the three DWD directions
in Figure \ref{figB-R:PanCanDWD}. A naive way to handle this is to
just plot the scores on the vertical and horizontal axes. An important
drawback to that is the scores plot no longer maintains relative positions
of the data in the underlying Euclidean space, and thus can give a
quite distorted view, especially when the direction vectors are far
from orthogonal. This is demonstrated by the toy data set in Figure
\ref{figF.H:NonOrthoAxes}. The data lie on four small circles in
$\mathbb{R}^{2}$, with more points on the circles centered on the
line from the blue to the cyan circles, shown using a solid line type.
Another direction of interest in such data is along the line from
the center of the red circle to the yellow one shown with a dot-dashed
line type. Recall that for both of those lines, projection of each
data point onto them consists of finding the closest point on the
line, which could be connected by a line segment that is perpendicular
to the given line, e.g the cyan line segments in the right panel of
Figure \ref{fig:2dToy-PC1} (these are not shown to keep the picture
from being too busy). The distance from the origin to each projection
is the corresponding score. Such scores are plotted in the other panels
of Figure \ref{figF.H:NonOrthoAxes}. The top right shows the red-yellow
scores (i.e. projections on the dot-dashed line) on the vertical axis,
and the blue-cyan scores on the horizontal axis. There is large distortion
in the relative position of the circles (the blue and red circles
appear much closer to each other, as do the cyan and yellow) and also
strong distortion in the shape of the circles. 

Such distortions generally occur from naively plotting scores with
respect to non-orthogonal axes\index{non-orthogonal axes}. Visualizations
that avoid such distortions are shown in the bottom two panels. The
bottom left has the same blue-cyan scores on the horizontal axis,
but the vertical axis is now the projections (i.e. scores) on the
line perpendicular to that axis, which is shown using the dashed line
type in the panels on the left. Everything is now a rigid rotation
of the quantities in the upper left panel, which thus preserves relative
locations. Such graphics were used in the plots above the diagonal
in Figure \ref{figB-R:PanCanDWD}. The lower right plot in Figure
\ref{figF.H:NonOrthoAxes} illustrates the plots below the diagonal.
There the horizontal axis represents the red-yellow direction (still
shown with the dot-dashed line type) and the vertical axis is the
orthogonal direction (again not shown as yet another line for ease
of understanding). As for the lower left panel there is no distortion
of the relative positions, and the projection of the other direction
of interest (the solid blue-cyan line) is shown to maintain the comparison
between the two directions. With this understanding in mind, it could
be useful to review the corresponding above and below diagonal plots
for the Pan Cancer data\index{Pan Cancer data} in Figure \ref{figB-R:PanCanDWD}.\index{scatterplot matrix|)}

\begin{figure}[H]
\includegraphics[scale=0.85]{ChapterF/OODAbookChpFFigHNonOrthoAxes}

\caption{Toy example showing problems with plotting projections onto non-orthogonal
lines. Raw 2-d data in upper left shows lines (directions of interest)
connecting clusters of interest. Upper right panel is a plot of projections
on red-yellow line (dot-dashed) as a function of projections on blue-cyan
line (solid), which introduces large distortions in placement of data.
Lower left panel shows projections on blue-cyan (solid) line, and
on orthogonal line (dashed), which keep relative positions, i.e. are
just rotations (and flips) of the original data. Lower right is same
for projections on the red-yellow (dot-dashed) line. \label{figF.H:NonOrthoAxes}}
\end{figure}
For comparison of slightly different populations, such as viewing
the impact of adding a few new observations, or partially perturbing
some data objects, a useful graphical device is to flip between two
different scatterplot matrices on a computer screen. This requires
use of exactly the same coordinate axes, which can be a challenge
for some implementations of PCA. The issue is that PCA is based on
an eigen analysis\index{eigen analysis|(} of the sample covariance
matrix (or equivalently on a Singular Value Decomposition\index{Singular Value Decomposition}
of the centered data matrix, as noted in Section \ref{subsecP-A-B:SVD}).
While eigen analysis provides a very useful decomposition, there is
a natural ambiguity in terms of the direction of the signs of eigenvectors
(direction vectors). In particular, for any eigenvectors all entries
can be multiplied by $-1$, i.e. the vector can be rotated $180^{\circ}$,
and it will still be a direction of maximal variation. This can strongly
disrupt attempts to flip between two graphics, and perhaps other tasks
as well, so it is useful to adopt a convention to minimize the occurrence
of this event. A recommended choice is to make all projections onto
the diagonal vector $\left(1,1,\cdots,1\right)^{t}$ (denoted $\boldsymbol{1}_{d,1}$
in (\ref{eq:ones matrix})) positive. This is simply achieved by choosing
the direction which makes the average of scores positive. This will
fail when the data have been row trait\index{row trait!mean|)}\index{trait|)}
mean centered, i.e. each column has mean 0 (since that inner product
gives the sum of the vector entries). In that case a reasonable alternative
is to choose the direction which makes the maximum absolute score
positive. The subspace generated by the diagonal vector $J_{d,1}$
is studied much more deeply in Section \ref{secP-A:PCA-Views}.\index{scores!scatterplot|)}

Note that many, but not all, of these scatterplot matrices are based
on PCA directions. Some exceptions are the use of DWD (e.g. the Pan
Can data in Figure \ref{figB-R:PanCanDWD}) to give directions with
better class separation. But there are many other potentially interesting
directions, depending on the context, as discussed in the next section.

\section{Alternatives to PCA Directions \label{sec F.D:AlternativeDirns}}

PCA\index{Principal Component Analysis!alternatives|(} is a well
proven workhorse in finding insightful directions for the computation
of scores that reveal structure in data, because the set of directions
of maximal variation which it targets frequently highlight such structure.
However PCA does not necessarily reveal all types of interesting structure,
especially when the latter appears in directions of relatively lower
variation. This issue becomes more serious in high dimensions. Hence
it is useful to realize that there are many alternatives, which can
give important insights depending on the particular context. These
include the following. Some of these come directly from multivariate
analysis (recall this is viewed as a special case of OODA at the beginning
of Chapter \ref{chap1A:What-is-OODA?}), but others arise from other
considerations.
\begin{itemize}
\item \emph{Classification Directions.\index{classification!directions}}
The statistical context called classification (also termed discrimination\index{discrimination})
is discussed in detail in Chapter \ref{chapJ:Classification}. Binary
linear classification methods seek to find a hyperplane that separates
two labeled subsets of the data usually in some optimal sense. The
normal vector to a separating plane often provides very useful directions
for visualization when there is interest in highlighting differences
between two subsets of data. The Mean Difference direction (computed
as the difference between the means of the two subsets, and discussed
in detail in Section \ref{secJ-A:ClassicalMethods}) is a simple and
generally interpretable direction for this purpose. Often better visual
impression comes from \index{Distance Weighted Discrimination|)}Distance
Weighted Discrimination (DWD studied in Section \ref{secJ-D:DWD}),
because the target of the optimization is more directly related to
visual separation than are subset means. Examples of this use of DWD
appear in Figures \ref{figB-R:PanCanDWD} and \ref{fig:HiDimToyDWD}.
In relatively low dimensional situations, good visual separation of
multiple classes can be achieved by \index{Canonical Variate Analysis}Canonical
Variate Analysis, see \citet{MardiaKentBibby1979MultivarAnal}, which
is also called Multiple Discriminant Analysis\index{Multiple Discriminant Analysis}
in \citet{DudaHartStork2001PattRecog}. Deeper discussion of related
issues can be found in Section \ref{secJ-D:DWD}.
\item \emph{Independent Component Analysis\index{Independent Component Analysis}.}
This source of interesting directions in high dimensional data space
seeks directions of maximal non-Gaussianity. The key idea was first
developed by \citet{Cardoso1993BlindSourceSep} and is discussed in
detail in \citet{Hyvarinen2001ICAbook}. Its main motivation is called
\emph{Blind Source Separation}\index{Blind Source Separation}, which
mathematically models the ability of human hearing to comprehend multiple
simultaneous conversations. This is done using two different mixtures
(linear combinations) of time series signals. When the coefficients
are known, signal separation only needs simple matrix inversion. More
challenging is the blind case of unknown coefficients, which Independent
Component Analysis solves by finding coefficients resulting in maximal
non-Gaussianity. This works because linear combination is essentially
an averaging operation that tends to increase Gaussianity. The approach
has become very popular in many neuroscience analyses. Such directions
are often useful in exploratory\index{exploratory analysis} data
analysis, as they tend to find interesting structure such as directions
of multi-modality. The idea of finding directions of interest by maximizing
non-Gaussianity goes back at least to the \emph{\index{varimax}varimax}
subspace rotation from \emph{factor analysis} proposed by \citet{Kaiser1958Varimax}.
\item \emph{\index{Fourier!subspace directions}Fourier Subspace.} For applications
where periodicity plays a central role, Fourier basis vectors can
provide useful directions for projection\index{projection|)}. This
idea was used to give an insightful analysis of the \citet{Spellman1998YeastCellCycle}
Yeast Cell Cycle data\index{Yeast Cell Cycle data} by \citet{ZhaoMarronWells2004YeastCellCycle}.
That noisy data set featured gene expression measured over two cell
cycles with a high level of background noise. The noise was high enough
that a standard PCA resulted in non-interpretable directions which
did not exhibit the periodicity that was the point of the study. A
much better noise reduction was done by projecting the data onto the
subspace of even Fourier frequencies, which kept only parts of the
data where the first half of the time series is identical to the second
half. Good access to aspects of Fourier analysis important to statistics
can be found in \citet{Brillinger1981} and \citet{Bloomfield2000}.
\item \emph{\index{Wavelet}Wavelet Decomposition.} This important variation
of Fourier Analysis shifts the focus from approximately periodic signals
to those with more local structure, resulting in orthogonal bases
with good signal compression (and thus efficient data object representation)
for curves with local structure such as jumps and bumps. See Section
\ref{secBA.2:OverviewObjFeat} for references and more discussion
of wavelets in the context of curve estimation. 
\item \emph{\index{known modes of variation}Known modes of variation\index{modes of variation!known}.}
In some applications particular modes of variation are given as very
important. For example in the \index{evolutionary biology|(}evolutionary
biology work of \citet{izemKingsolver2005Variation}, the focus is
on modes of variation of central importance to the ability of a population
to adapt to environmental changes. Some of these modes were nonlinear,
motivating some relatively early research on quantifying variation
on manifold\index{manifold}s. Theoretical properties of the methods
developed there were investigated in \citet{IzemMarron2007NonLinearVariation}.
More recent developments in the area on nonlinear modes of variation
are discussed in Chapters \ref{chapG:Manifold-Data-Analysis} and
\ref{chapH:FDA-Curve-Registration}.
\item \emph{\index{Maximal Smoothness Directions}Maximal Smoothness Directions.}
As noted several times above PCA is very powerful because of, but
can also sometimes be hobbled by, its property of finding directions
of maximal variation in a Euclidean feature space\emph{\index{feature space}}.
In some evolutionary biology\index{evolutionary biology|)} applications,
as discussed in \citet{GaydosKingsolver2013SmoothDirections} and
\citet{Kingsolver2015EvolConstraints}, the major interest is instead
on directions of \emph{minimum variation}. This is because they represent
ways in which a given population is \emph{unable} to adapt to environmental
change, called the \emph{\index{nearly null space}nearly null space}.
Such directions can also be easily found by PCA, simply by working
with small eigenvalues, instead of the large ones that are the basis
of conventional PCA. However, unlike PCA directions, which are typically
interpretable because they are driven by important modes of variation,
directions of minimal variation are generally far less interpretable
as most of them are driven by uninteresting sampling artifacts. To
find important interpretable directions of minimal variation, among
the nearly null subspace that minimizes variation, those papers propose
finding directions of \emph{maximal smoothness}. These are found by
eigen analysis\index{eigen analysis|)} of a matrix reflecting \emph{simplicity}
(essentially smoothness, but also usually reflecting maximal interpretability)
of directions instead of variation as is done in conventional PCA.
In particular, \citet{GaydosKingsolver2013SmoothDirections} propose
a method called \emph{PrinSimp.} This studies PCA eigenvectors from
a perspective that simultaneously accounts for variation explained
as well as simplicity of representation of the null space. \citet{Kingsolver2015EvolConstraints}
provide a number of interesting applications of \index{PrinSimp method}PrinSimp.
An R software package implementation of PrinSimp is available from
\citet{ZhangHeckman2014prinsimp}.
\item \emph{Non-Negative Matrix Factorization\index{Non-Negative Matrix Factorization}.}
In some functional data applications, not only are the input data
curves nonnegative, but there is also a desire for all aspects of
the analysis to be nonnegative. An important example is chemical spectra
as data objects, where the height of the curves represent mass, for
example the TIC curves in Figure \ref{fig:TIC-spectra}. Standard
functional PCA can be very unattractive to chemists, because while
the first component direction will typically have nonnegative loadings,
the orthogonality constraint of PCA entails that all other components
will have negative loadings, which chemists tend to find challenging
to interpret (because mass is an intrinsically nonnegative quantity).
Such considerations are addressed by \emph{Non-negative Matrix Factorization}
(NMF) as proposed by \citet{LeeSeung1999NNMF}. Simple insight into
the relationship between NMF and PCA comes from considering (as detailed
in Sections \ref{chap1A:What-is-OODA?} and \ref{subsecP-A-B:SVD})
the latter as providing an optimal approximation of a mean centered\index{mean!center|)}
$d\times n$ data matrix $\boldsymbol{X}$, by a matrix product of
a $d\times r$ loadings matrix $\boldsymbol{L}$ (e.g. in the Twin
Arches\index{Twin Arches data|)} analysis of Figure \ref{fig:50d-ToyPCA}
the lower left modes of variation are multiples of the $r=3$ columns
of $\boldsymbol{L}$) and an $r\times n$ scores matrix $\boldsymbol{S}$
(e.g. in Figure \ref{fig:50d-ToyPCAScatMat} various pairs of the
$r=4$ rows of $\boldsymbol{S}$ are plotted against each other).
In particular, as seen in Section \ref{secP-A:PCA-Views} $\boldsymbol{L}$
and $\boldsymbol{S}$ essentially provide minimizers (subject to orthogonality
conditions) of the Frobenius norm\index{Frobenius norm} (defined
at (\ref{eq:Frobeniusnorm})) of $\boldsymbol{X}-\boldsymbol{L}\boldsymbol{S}$.
NMF is based on a similar optimization, but with the constraints that
all the entries of both $\boldsymbol{L}$ and $\boldsymbol{S}$ are
non-negative. As noted say in \citet{FogelYoung2006NMFmicroarrays}
there are other applications, e.g. in gene expression, where an essentially
additive, nonnegative decomposition of variation can be viewed as
more intuitive than the orthogonal decomposition provided by PCA.
As discussed in Section \ref{secG-C:BackwardsPCA}, a weakness of
classical NMF is that it is not \emph{nested}\index{nested analyses}.
This issue has been addressed by the Nested Nonnegative Cone Analysis
proposed by \citet{ZhangLuMarron2015NNCA}. 
\item \emph{Binary Matrix Factorization\index{Binary Matrix Factorization}.}
In some applications, the data matrix $\boldsymbol{X}$ may be \emph{binary},
i.e. each entry may take on only the values $0$ or $1$. The entries
are numbers, so one can analyze the variation using a classical PCA\index{Principal Component Analysis!alternatives|)}
decomposition. But for most purposes, such a decomposition will not
be very interpretable, because neither loadings nor scores\index{scores|)}
will be binary. A more interpretable version is offered by \emph{Binary
Matrix Factorization} proposed by \citet{Srebro2005BMF} and \citet{Zhang2007bBMF},
which has a goal in the spirit of NMF, i.e. an approximation of the
form $d\left(\boldsymbol{X},\boldsymbol{L}\boldsymbol{S}\right)$,
where the matrix product $\boldsymbol{L}\boldsymbol{S}$ is computed
with Boolean operations\index{Boolean operations}, and $d$ is an
appropriate metric\index{metric} on binary matrices such as \citet{Hamming1950ErrorCorrectingCodes}
distance which just counts the number of different entries. Also appealing
is the \citet{Jaccard1901Index} distance which is the proportion
of the number of differences among the total number of ones.
\item \emph{Lines Containing Data.} In computer vision, when analyzing images
involving man-made objects, such as touristic city scapes, many useful
methods are based on finding lines in an image. This is done by filtering
the image to find \emph{points of interest}, e.g. points of maximal
change in various senses. Lines are then discovered by finding lines
containing a reasonable number of such points, which are used as building
blocks for more complex types of analysis. The RANSAC (more precisely
Random Sample Consensus\index{Random Sample Consensus}) method, proposed
by \citet{FischlerBolles1981RANSAC}, is a useful algorithm for finding
such \emph{inlier lines}\index{inlier lines}, featuring a few but
``very close to collinear'' points. These ideas are an interesting
inversion of the robust\index{robust methods}ness ideas of Chapter
\ref{chapO:Robust-Methods} where the focus is on outliers.
\item \emph{\index{Maximal Data Piling}Maximal Data Piling.} A perhaps
surprising phenomenon in high dimensional data analysis (when the
dimension is larger than the sample size), is that given two subsets
of the data, there is (with high probability) a set of directions
where each subset projects to a single point. The direction which
maximizes the distance between these points is called the \emph{Maximal
Data Piling} direction, which is studied in more detail in Section
\ref{secJ-A:ClassicalMethods}. This direction is not particularly
useful for most visualization purposes, but is still worth keeping
in mind when contemplating the large variety of directions available
in high dimensional space. 
\item \emph{Auxiliary Data Directions.} In many situations, visualization
in directions that incorporate additional data are very useful. The
simplest of these may be standard multiple linear regression, when
the additional information is in terms of a response value $Y_{i}$
for each data object. More sophisticated approaches, especially in
high dimensions include the \emph{Sliced Inverse Regression\index{Sliced Inverse Regression}}
of \citet{Li1991SIR} and the \emph{Sufficient Dimension Reduction\index{Sufficient Dimension Reduction}}
of  \citet{CookLee1999SDR}. When richer vector-valued data are available,
methods such as \index{Partial Least Squares}Partial Least Squares
and \index{Canonical Correlation Analysis}Canonical Correlation Analysis,
discussed in Section \ref{secP-B:TwoBlockDecomp}, also provide very
useful directions for visualization\index{visualization|)}.\index{modes of variation|)}
\end{itemize}

