在該方法中影響聚類結果的主要因素包括資料預處理方式、距離測度和聚類演算法的選擇。
聚類演算法表示一組決定將樣本分配給哪個組（聚類群）的規則。可以為該過程打個比方，在一個滿是人（樣本）的酒吧中，決定坐在哪一張桌子上（聚類群、組），有些人是我的好朋友（低距離，或高相似性），我傾向和他們坐在一起；有些人我比較陌生（高距離，或低相似性），我通常和他們離的較遠。
接下來簡介層次聚合分類及其在R語言中的計算方法。
層次聚合分類方法概述
首先簡介幾種常見的層次聚合分類模式。
單連線聚合聚類（單聯動）
單連線聚合聚類（single linkage agglomerative clustering），或稱最近鄰體分類（nearest neighbour sorting），依據最短的成對距離（或最大的相似性）聚合物件（Lukaszewicz, 1951）。
以“坐在酒吧裡的每個人”為例，都更傾向於和最近的朋友們坐在一桌（可理解為相似度最高或距離最短），而不考慮和陌生人湊在一桌（可理解為相似度較低或距離較大）。
首先計算相似度矩陣，然後根據物件對的組成相似度對物件進行排序。最相似的對將形成第一個簇，即樣本212與樣本214最相似歸入第一簇，在S20=0.600處合併；排第二位的最相似對將形成第二個簇，即樣本431與樣本432歸入第二簇，在S20=0.500處合併；以此類推，直到將所有物件加入簇中。
因此在單連線聚合聚類中，兩個聚類群之間的距離定義為每個簇中兩個物件之間的最短距離，這種計算規則使聚類很容易實現。
不難看出，單連線聚合聚類使一個物件很容易聚合到一個分類組中，因為一次連線足以導致融合，因此單連線聚合聚類也被稱為“最親密朋友”（closest friend）法。並且單連線聚類產生的梯度非常明顯，即容易識別資料的梯度。
單連線聚類由於其是空間壓縮的，這種演算法使相異較大的物件很快合併在一起。因此其缺點是產生的分類組之間的差異通常不清晰，即從分割槽的角度來看比較難以解讀，並且錯分的可能性也較大。
完全連線聚合聚類（全聯動）
與單連線聚合聚類相反，完全連線聚合聚類（Complete linkage agglomerative clustering），或稱最遠鄰體分類（Furthest neighbour sorting），允許一個物件（或一個組）與另一個組聚合的依據是最遠距離對，即一個組吸納一個新成員的依據是比較該組內所有成員與備選新成員的最遠距離（Lance and Williams, 1967）。在該模式下，兩個組所有成員之間的距離都必須全部計算，然後再比較。
繼續以“坐在酒吧裡的每個人”舉例，完全連線聚類可以理解為我加入了一個具有我不喜歡的人的桌子中，但是在我加入後該桌的關係也並不糟糕，因為其它每個桌子上都有我更不喜歡的人。
完全連線聚合聚類中，兩個聚類群之間的距離定義為每個簇中兩個物件之間的最長距離。
在完全連線聚類中，由於合併後的組與其它組的距離是原來距離的最大者，加大了並組與其它組間的距離，因此其是空間擴張的，對進一步合併起著阻礙作用。對於已知組，吸收一個新成員要求該組內所有成員都參與比較，組越大就越難再加入新成員。
因此與單連線聚類方法相比，完全連線聚類更傾向產生很多小的分離的組，並且分類組之間的差異也比較明顯，其更適合尋找和識別資料的間斷分佈。
平均聚合聚類
單連線聚合聚類或完全連線聚合聚類均屬於基於連線的聚類方法。
與上述方法相比，平均聚合聚類（average agglomerative clustering）是一類基於物件間平均相異性或聚類簇形心（centroid）的聚類方法，目的將物件加入到其平均距離最低的組。在“坐在酒吧裡的每個人”例子中，平均聚合聚類模式可以理解為我綜合衡量所有的桌中，對每張桌子上的人的平均好感度，然後加入具有最高平均好感度的那桌。
平均聚合聚類中，兩個聚類群之間的距離定義為一個簇中每個物件到另一個簇中每個物件之間的平均距離。
平均聚合聚類包含4種主要型別，區別在於組的位置計算方式（中值或質心）和當計算融合距離時是否用每組包含的物件數量作為權重（Sneath and Sokal, 1973）。
下圖形象展示了單連線聚類（a）、完全連線聚類（b）、中值聚類（c）和質心聚類（d）的區別（Greig-Smith 1983）。
對於中值法，一個物件加入一個組的依據是這個物件與該組每個成員之間的平均距離，兩個組聚合的依據是一個組內所有成員與另一組內成員之間所有物件對的平均距離。中值法保持了合併後的組與其它物件（或組）間的距離，即空間保持的，因此通常比基於連線的聚類方法（單連線聚類或完全連線聚類）效果更好；缺點是非單調性（Gauch, 1982）。
質心法將兩個組間的距離定義為兩個組質心間的距離，組質心指組內所有物件間距離的均值。在單個物件或只有兩個物件組成的組中，質心法與中值法完全相同，但當物件數大於3時，二者將不同。質心法的缺點在於，當兩個組所含物件數目相差懸殊時，新合併的組將更接近較大的一組，使較小組的特徵消失，而影響聚類結果（Greig-Smith, 1983）。
此外質心法還可能會導致聚類樹翻轉的問題，詳見下文R語言執行部分。
Ward最小方差聚類
Ward最小方差聚類（Ward’s minimum variance clustering）是一種基於最小二乘法線性模型準則的聚類方法，分組的依據是使組內誤差平方和（每個物件到其組質心的距離的平方和）最小化（Ward,1963）。
Ward聚類方法試圖尋找一種解決方案，將兩組物件（簇）合併在一起。一開始，每個物件形成一個獨立的簇，誤差平方和為零。在隨後的每一步中，該方法將物件與物件、物件與簇或簇與簇進行合併，此時誤差平方和增大，並最終找到具有最小誤差平方和的組合，使總體誤差平方和的增幅最低。
圖示兩種不同的Ward聚類計算過程。對於每個聚類簇k，計算單個物件（yij(k)）到該簇質心（mj(k)）的距離的平方和（左式），或者取簇內所有物件對之間的平均平方距離（Dhi2，右式）。
由於物件到質心的距離是在歐幾里得空間中計算的，因此在該方法中只能使用歐式距離測度（如歐幾里得距離）。對於其它非歐式屬性的距離測度的使用，可能會使聚類結果違背最小化組內方差的基本原則，使聚類結果不可靠。儘管如此，可以透過轉化的方式實現它們的應用，例如Bray-curtis距離，可在平方根轉化後獲得歐式幾何特徵。
靈活聚類
Lance和Williams（1966，1967）提出了一種能夠涵蓋上述聚類方法的模型，稱為靈活聚類（flexible clustering)，該方法中沒有固定的聚合策略，透過指定不同的引數項（αh、αi、β、γ）實現不同的聚類效果。
式中，g、h、i代表了三個聚類群（組），S為組間相似性，n為各聚類群內的物件數量。
若使用距離測度而非相似性，則可使用下式計算，D代表了組間距離。
對於上述所提到的幾種聚類方法，在靈活聚類中均可視為特殊形式。
可知，其中最主要的值是β，因此大多數情況下，靈活聚類也稱為靈活beta連線聚類（flexible beta linkage clustering）。
β稱為聚集強度係數（cluster-intersing coefficient），取值範圍[-1,1)，當β=-1時，具有最低水平的聚類程度，各物件單獨成簇；當β越接近於1，越具更高水平的聚類程度，或者說更多的物件傾向於聚為一簇。
圖示基於Bray-curtis距離的靈活聚類，3種不同的結果分別對應於β=-0.95（左）、β=-0.25（中）和β=0.95（右）時的情形。
層次聚合分類在R語言中的計算
接下來展示上述聚類方法在R語言中的實現。
例如期望透過聚類識別不同組織的基因表達相似性，不同環境群落物種組成相似性等。
下文的示例資料，R程式碼等可見百度盤：
https://pan.baidu.com/s/1tM2iSX6-r0aWxZmkXD89cA
計算距離
計算樣本（或分組）間的關聯矩陣（或距離矩陣）是進行層次聚合分類分析的第一步，因此選擇恰當的關聯絡數非常重要。
##示例 1，基因表達資料
dat_gene <- read.delim('gene_express.txt', row.names = 1)
dat_gene <- t(dat_gene)
#常用歐幾里得距離表示各樣本中的基因表達相異性，可由 vegan 包 vegdist() 計算距離測度
dis_euc <- vegan::vegdist(dat_gene, method = 'euclidean')
##示例 2，物種丰度資料
dat_otu <- read.delim('otu_table.txt', row.names = 1)
dat_otu <- t(dat_otu)
#群落資料的分析中，Bray-curtis 距離是常用的距離測度型別
dis_bray <- vegan::vegdist(dat_otu, method = 'bray')
#不過考慮到某些聚類方法只能以歐式距離作為輸入
#因此可以透過平方根轉化的形式，將 Bray-curtis 距離轉化為具有歐式幾何屬性的距離測度
ade4::is.euclid(dis_bray)
dis_bray_euc <- dis_bray^0.5
ade4::is.euclid(dis_bray_euc)
##此外，如果提供了現有的已經計算好的某種距離矩陣，也可直接讀取
#例如，外部的 Bray-curtis 距離矩陣
dis_bray <- as.dist(read.delim('bray_distance.txt', row.names = 1))
#平方根轉化，使其具有歐式幾何屬性
dis_bray_euc <- dis_bray^0.5
執行層次聚合分類
單連線聚合聚類&完全連線聚合聚類
可使用hclust()計算。
##基於連線的層次聚合分類，詳情 ?hclust
#以物種丰度的 Bray-curtis 距離為例
#單連線聚合聚類，單聯動
clust_single <- hclust(dis_bray, method = 'single')
summary(clust_single)
#完全連線聚合聚類，全聯動
clust_complete <- hclust(dis_bray, method = 'complete')
summary(clust_complete)
#聚類樹，預設效果
par(mfrow = c(1, 2))
plot(clust_single, main = 'single')
plot(clust_complete, main = 'complete')
平均聚合聚類
可使用hclust()計算。
##平均聚合聚類，詳情 ?hclust
#以物種丰度的 Bray-curtis 距離為例
#使用中值的非權重成對組法，UPGMA
clust_average <- hclust(dis_bray, method = 'average')
summary(clust_average)
#使用質心的非權重成對組法，UPGMC
clust_centroid <- hclust(dis_bray, method = 'centroid')
summary(clust_centroid)
#使用中值的權重成對組法，WPGMA
clust_mcquitty <- hclust(dis_bray, method = 'mcquitty')
summary(clust_mcquitty)
#使用質心的權重成對組法，WPGMC
clust_median <- hclust(dis_bray, method = 'median')
summary(clust_median)
#聚類樹，預設效果
par(mfrow = c(2, 2))
plot(clust_average, main = 'UPGMA')
plot(clust_centroid, main = 'UPGMC')
plot(clust_mcquitty, main = 'WPGMA')
plot(clust_median, main = 'WPGMC')
質心法的聚類方法有時會導致聚類樹翻轉（reversal），如上UPGMC或WPGMC所示，聚類樹不再形成連續的巢狀分割槽，使分類結果難以解讀。
對於該現象產生的成因，可簡要描述如下。存在3個物件X1、X2和X3（a），X1和X2首先歸為一簇，因為它們之間具有最短的距離；儘管X3-X1或X3-X2的距離大於X1-X2的距離，但是X3到X1和X2質心的距離小於X1-X2的距離，此時即產生聚類樹翻轉的情形（b）。
Legendre和Legendre（1998）在“Numerical Ecology”第342頁前後對該現象作了詳細描述，並建議用多分法（polychotomy）代替傳統的二分法（dichotomy）解讀這種圖。
Ward最小方差聚類
可使用hclust()計算。
##Ward 最小方差聚類（包含兩種不同的 Ward 聚類演算法，詳情 ?hclust）
#該方法中只能使用歐式幾何屬性的距離測度，因此更換為平方根轉化後的 Bray-curtis 距離測度
#ward.D
clust_ward <- hclust(dis_bray_euc, method = 'ward.D')
summary(clust_ward)
#ward.D2
clust_ward2 <- hclust(dis_bray_euc, method = 'ward.D2')
summary(clust_ward2)
#聚類樹，預設效果
par(mfrow = c(1, 2))
plot(clust_ward, main = 'ward.D')
plot(clust_ward2, main = 'ward.D2')
靈活聚類
可使用cluster包的agnes()計算。
##靈活聚類，詳情 ?agnes
#也可在該函式中透過指定 method 執行上述的多種常規聚類，如 method='average' 即為 UPGMA
#當 method='flexible' 時，par.method 引數中的 4 個值分別對應了 αh、αi、γ、β 的設定值，例如
clust_flex <- cluster::agnes(dis_bray_euc, method = 'flexible', par.method = c(0.5, 0.5, 0, 0.5))
summary(clust_flex)
plot(clust_flex)
聚類評估
層次聚類是探索性分析，而非統計檢驗的過程。
對於如何比較和評估層次聚類結果，尋找可解讀的聚類簇等，放在下篇介紹。
聚類樹視覺化調整
以下再簡單展示一些聚類樹的調整示例。
同樣地，聚類樹更詳細的視覺化方法將在後文使用一篇單獨展示。
#plot() 繪製聚類樹引數簡要說明，以上述 UPGMA 結果為例
#預設
plot(clust_average)
#hang 引數（預設 0.1）可用於調整標籤與樹高的關係，當hang < 0時，強制標籤從 0 開始並對齊，如下示例
par(mfrow = c(1, 2))
plot(clust_average, hang = 0.2)
plot(clust_average, hang = -1)
#main、sub、xlab、ylab 用於修改圖片標題或座標軸標題，如下示例
plot(clust_average, main = 'UPGMA\n(Bray-curtis distance)', sub = '', xlab = 'Sample', ylab = 'Height')
#col 可定義顏色，cex 可調整標籤字型大小，lwd 可調整線寬，lty 可調整線型，如下示例
plot(clust_average, col = 'red', cex = 1, lwd = 1, lty = 2)
#設定選擇觀測的聚類群數量為 4
plot(clust_average, main = 'UPGMA\nBray-curtis distance', sub = '', xlab = 'Sample', ylab = 'Height')
rect.hclust(clust_average, k = 4, border = c('red', 'blue', 'green3', 'orange'))
#vegan 包函式 reorder() 可以實現對聚類樹重排，詳情 ?reorder
#儘可能使聚類樹內物件的排列順序和原始相異矩陣內物件的排列順序一致
#不過效果似乎一般
library(vegan)
plot(reorder(clust_average, dis_bray), main = 'UPGMA\nBray-curtis distance', sub = 'after reorder', xlab = 'Sample', ylab = 'Height')
rect.hclust(reorder(clust_average, dis_bray), k = 4, border = c('red', 'blue', 'green3', 'orange'))
參考資料
