I am performing hierarchical clustering on data I've gathered and processed from the reddit data dump on Google BigQuery. 

My process is the following:


Get the latest 1000 posts in /r/politics
Gather all the comments
Process the data and compute an n x m data matrix (n:users/samples, m:posts/features)
Calculate the distance matrix for hierarchical clustering
Choose a linkage method and perform the hierarchical clustering
Plot the data as a dendrogram


My question is, how do I determine what the best linkage method is? I'm currently using Ward but how do I know if I should be using single, complete, average, etc?

I'm very new to this stuff but I can't find a clear answer online as I'm not sure there is one. So what might be a good idea for my application? Note the data is relatively sparse in the sense that the n x m matrix has a lot of zeroes (most people don't comment on more than a few posts).

Methods overview
Short reference about some linkage methods of hierarchical agglomerative cluster analysis (HAC).
Basic version of HAC algorithm is one generic; it amounts to updating, at each step, by the formula known as Lance-Williams formula, the proximities between the emergent (merged of two) cluster and all the other clusters (including singleton objects) existing so far. There exist implementations not using Lance-Williams formula. But using it is convenient: it lets one code various linkage methods by the same template.
The recurrence formula includes several parameters (alpha, beta, gamma). Depending on the linkage method, the parameters are set differently and so the unwrapped formula obtains a specific view. Many texts on HAC show the formula, its method-specific views and explain the methods. I would recommend articles by Janos Podani as very thorough.
The room and need for the different methods arise from the fact that a proximity (distance or similarity) between two clusters or between a cluster and a singleton object could be formulated in many various ways. HAC merges at each step two most close clusters or points, but how to compute the aforesaid proximity in the face that the input proximity matrix was defined between singleton objects only, is the problem to formulate.
So, the methods differ in respect to how they define proximity between any two clusters at every step. "Colligation coefficient" (output in agglomeration schedule/history and forming the "Y" axis on a dendrogram) is just the proximity between the two clusters merged at a given step.

Method of single linkage or nearest neighbour. Proximity
between two clusters is the proximity between their two closest
objects. This value is one of values of the input matrix. The conceptual metaphor of this build of cluster, its archetype, is spectrum or chain. Chains could be straight or curvilinear, or could be like "snowflake" or "amoeba" view. Two most dissimilar cluster members can happen to be very much dissimilar in comparison to two most similar. Single linkage method controls only nearest neighbours similarity.

Method of complete linkage or farthest neighbour. Proximity
between two clusters is the proximity between their two most distant
objects. This value is one of values of the input matrix. The metaphor of this build of cluster is circle (in the sense, by hobby or plot) where two most distant from each other members cannot be much more dissimilar than other quite dissimilar pairs (as in circle). Such clusters are "compact" contours by their borders, but they are not necessarily compact inside.

Method of between-group average linkage (UPGMA). Proximity
between two clusters is the arithmetic mean of all the proximities
between the objects of one, on one side, and the objects of the
other, on the other side. The metaphor of this build of cluster is quite generic, just united class or close-knit collective; and the method is frequently set the default one in hierarhical clustering packages. Clusters of miscellaneous shapes and outlines can be produced.

Simple average, or method of equilibrious between-group average linkage (WPGMA) is the modified previous. Proximity between two clusters is the
arithmetic mean of all the proximities between the objects of one, on
one side, and the objects of the other, on the other side; while the
subclusters of which each of these two clusters were merged recently
have equalized influence on that proximity – even if the subclusters
differed in the number of objects.

Method of within-group average linkage (MNDIS). Proximity between
two clusters is the arithmetic mean of all the proximities in their
joint cluster. This method is an alternative to UPGMA. It usually will lose to it in terms of cluster density, but sometimes will uncover cluster shapes which UPGMA will not.

Centroid method (UPGMC). Proximity between two clusters is the proximity between their geometric centroids: [squared] euclidean
distance between those. The metaphor of this build of cluster is proximity of platforms (politics). Like in political parties, such clusters can have fractions or "factions", but unless their central figures are apart from each other the union is consistent. Clusters can be various by outline.

Median, or equilibrious centroid method (WPGMC) is the modified previous. Proximity between two clusters is the proximity between their geometric
centroids ([squared] euclidean distance between those); while the
centroids are defined so that the subclusters of which each of these
two clusters were merged recently have equalized influence on its
centroid – even if the subclusters differed in the number of objects. Name "median" is partly misleading because the method doesn't use medians of data distributions, it is still based on centroids (the means).

Ward’s method, or minimal increase of sum-of-squares (MISSQ), sometimes incorrectly called "minimum variance" method. Proximity
between two clusters is the magnitude by which the summed square in
their joint cluster will be greater than the combined summed square
in these two clusters: $SS_{12}-(SS_1+SS_2)$. (Between two singleton objects
this quantity = squared euclidean distance / $2$.) The metaphor of this build of cluster is type. Intuitively, a type is a cloud more dense and more concentric towards its middle, whereas marginal points are few and could be scattered relatively freely.


Some among less well-known methods (see Podany J. New combinatorial clustering methods // Vegetatio, 1989, 81: 61-77.) [also implemented by me as a SPSS macro found on my web-page]:

Method of minimal sum-of-squares (MNSSQ). Proximity between two
clusters is the summed square in their joint cluster: $SS_{12}$. (Between
two singleton objects this quantity = squared euclidean distance /
$2$.)

Method of minimal increase of variance (MIVAR). Proximity between
two clusters is the magnitude by which the mean square in their joint
cluster will be greater than the weightedly (by the number of
objects) averaged mean square in these two clusters:
$MS_{12}-(n_1MS_1+n_2MS_2)/(n_1+n_2) = [SS_{12}-(SS_1+SS_2)]/(n_1+n_2)$. (Between two
singleton objects this quantity = squared euclidean distance / $4$.)

Method of minimal variance (MNVAR). Proximity between two
clusters is the mean square in their joint cluster: $MS_{12} =
   SS_{12}/(n_1+n_2)$. (Between two singleton objects this quantity = squared
euclidean distance / $4$.).


Still other methods represent some specialized set distances. HAC algorithm can be based on them, only not on the generic Lance-Williams formula; such distances include, among other: Hausdorff distance and Point-centroid cross-distance (I've implemented a HAC program for SPSS based on those.)
First 5 methods described permit any proximity measures (any similarities or distances) and results will, naturally, depend on the measure chosen.
Next 6 methods described require distances; and fully correct will be to use only squared euclidean distances with them, because these methods compute centroids in euclidean space. Therefore distances should be euclidean for the sake of geometric correctness (these 6 methods are called together geometric linkage methods). At worst case, you might input other metric distances at admitting more heuristic, less rigorous analysis. Now about that "squared". Computation of centroids and deviations from them are most convenient mathematically/programmically to perform on squared distances, that's why HAC packages usually require to input and are tuned to process the squared ones. However, there exist implementations - fully equivalent yet a bit slower - based on nonsquared distances input and requiring those; see for example "Ward-2" implementation for Ward's method. You should consult with the documentation of you clustering program to know which - squared or not - distances it expects at input to a "geometric method" in order to do it right.
Methods MNDIS, MNSSQ, and MNVAR require on steps, in addition to just update the Lance-Williams formula, to store a within-cluster statistic (which depends on the method).
Methods which are most frequently used in studies where clusters are expected to be solid more or less round clouds, - are methods of average linkage, complete linkage method, and Ward's method.
Ward's method is the closest, by it properties and efficiency, to K-means clustering; they share the same objective function - minimization of the pooled within-cluster SS "in the end". Of course, K-means (being iterative and if provided with decent initial centroids) is usually a better minimizer of it than Ward. However, Ward seems to me a bit more accurate than K-means in uncovering clusters of uneven physical sizes (variances) or clusters thrown about space very irregularly. MIVAR method is weird to me, I can't imagine when it could be recommended, it doesn't produce dense enough clusters.
Methods centroid, median, minimal increase of variance – may give sometimes the so-called reversals: a phenomenon when the two clusters being merged at some step appear closer to each other than pairs of clusters merged earlier. That is because these methods do not belong to the so called ultrametric. This situation is inconvenient but is theoretically OK.
Methods of single linkage and centroid belong to so called space contracting, or “chaining”. That means - roughly speaking - that they tend to attach objects one by one to clusters, and so they demonstrate relatively smooth growth of curve “% of clustered objects”. On the contrary, methods of complete linkage, Ward’s, sum-of-squares, increase of variance, and variance commonly get considerable share of objects clustered even on early steps, and then proceed merging yet those – therefore their curve “% of clustered objects” is steep from the first steps. These methods are called space dilating. Other methods fall in-between.
Flexible versions. By adding the additional parameter into the Lance-Willians formula it is possible to make a method become specifically self-tuning on its steps. The parameter brings in correction for the being computed between-cluster proximity, which depends on the size (amount of de-compactness) of the clusters. The meaning of the parameter is that it makes the method of agglomeration more space dilating or space contracting than the standard method is doomed to be. Most well-known implementation of the flexibility so far is to average linkage methods UPGMA and WPGMA (Belbin, L. et al. A Comparison of Two Approaches to Beta-Flexible Clustering // Multivariate Behavioral Research, 1992, 27, 417–433.).
Dendrogram. On a dendrogram "Y" axis, typically displayed is the proximity between the merging clusters - as was defined by methods above. Therefore, for example, in centroid method the squared distance is typically gauged (ultimately, it depends on the package and it options) - some researchers are not aware of that. Also, by tradition, with methods based on increment of nondensity, such as Ward’s, usually shown on the dendrogram is cumulative value - it is sooner for convenience reasons than theoretical ones. Thus, (in many packages) the plotted coefficient in Ward’s method represents the overall, across all clusters, within-cluster sum-of-squares observed at the moment of a given step. Don't miss to read the documentation of your package to find out in which form the particular program displays colligation coefficient (cluster distance) on its dendrogram.
One should refrain from judging which linkage method is "better" for his data by comparing the looks of the dendrograms: not only because the looks change when you change what modification of the coefficient you plot there - as it was just described, - but because the look will differ even on the data with no clusters.
To choose the "right" method
There is no single criterion. Some guidelines how to go about selecting a method of cluster analysis (including a linkage method in HAC as a particular case) are outlined in this answer and the whole thread therein.
