生長曲線構建的不同方法比較*

劉 裕何健榮郭 勇夏曉燕王 平莫偉健盧錦華李偉棟于 佳林穗方馮 瓊邱 琇△

生長曲線構建的不同方法比較*

劉 裕1，2何健榮1郭 勇1夏曉燕1王 平1莫偉健1盧錦華1李偉棟1于 佳1林穗方1馮 瓊1邱 琇1△

目的尋找適合于構建人體測量指標生長標準曲線的方法和途徑。方法以構建廣州市胎齡別新生兒出生體重百分位曲線為例，從Tukey方法、穩健性回歸和高斯混合模型等三種異常數據識別和剔除方法中獲得最佳數據預處理效果，然后對比三次樣條方法、LMS方法和GAMLSS方法對百分位數曲線構建的影響。結果高斯混合模型對多峰分布數據中的主要分布識別比較理想，而對單峰分布，穩健性回歸比Tukey方法更加可靠。而從擬合優度以及小于胎齡兒（SGA）、大于胎齡兒（LGA）的識別能力看，GAMLSS構建的胎齡別新生兒出生體重百分位曲線比三次樣條和LMS方法估計精度更高。結論數據預處理過程應根據數據分布的特點選用合適的異常值識別和剔除方法，而曲線光滑過程中，GAMLSS方法可以對四階矩進行建模，得到的百分位數曲線平滑且誤差更小。

生長曲線 Tukey方法 穩健性回歸 高斯混合模型 LMS方法 GAMLSS

醫學參考值范圍是臨床實踐中篩查異常情況的重要依據。當測量指標依賴于其他協變量（如年齡）時，應制定隨協變量變化的參考值曲線，例如兒童青少年的身高、體重、體質指數、血壓、肺功能等生長參考值曲線［1－4］。構建參考值曲線的關鍵是對原始數據進行預處理并采用合適的統計學方法描繪曲線。盡管目前有多種統計學方法可用于曲線的構建，但不同的數據處理和統計學方法對曲線構建產生的影響，鮮有文獻報道。本文擬利用廣州市新生兒出生體重數據，對此進行探討，并提出方法選擇的建議。

資料來源

本研究的數據來源于廣州市圍產保健－產時管理軟件。該系統自2000年啟用，由廣州所有助產機構錄入在院內分娩產婦及新生兒的出生及圍產期信息，并用于核發出生醫學證明。數據質量控制每季度由市區兩級婦幼保健院完成。本研究納入2009年1月至2011年12月出生胎齡在26～43周的單胎活產男孩資料。考慮到26～33周胎齡兒數量較少，將2007年1月至2008年12月期間出生的26～33周胎齡兒納入分析。

方 法

1．胎齡的確定

依據廣州市助產技術管理規范，新生兒胎齡用妊娠16周前B超檢查測定的頭臀徑或雙頂徑進行確定；無法獲得B超檢查結果者，胎齡由末次月經結合出生后胎齡評估確定。胎齡以整周表示。

2．異常值識別

胎齡的錯分和出生體重的極端值會給胎齡別出生體重百分位數曲線的構建帶來偏移［3－7］，具體表現為，胎齡別出生體重呈不對稱、雙峰或長尾分布，需要首先識別和剔除。針對不同的分布特征，常用的識別異常數據的穩健性統計方法包括Tukey法［5－6］、穩健性回歸（Robust regression）［3－4］和高斯混合模型（Gaussian m ixturemodel）［3，7］。為了比較不同異常數據處理方法對生長曲線構建的影響，分別采用Tukey法和高斯混合模型進行數據預處理，計算胎齡別出生體重的第10、50和90百分位數曲線。由于高斯混合模型僅用于多峰分布的情形，對35周后的體重數據采用穩健性回歸方法。三種方法的原理如下：

Tukey方法 對每個胎齡x的出生體重數據y，首先計算第1四分位數（Q1）、第3四分位數（Q3）及其四分位數間距（IQR＝Q3－Q1）；然后計算L統計量，以反映兩端的體重數值距離Q3或Q1的距離是IQR的多少倍：（1）若y＞Q3，L＝（y－Q3）／IQR；（2）若y∈［Q1，Q3］，L＝0；（3）若y＜Q1，L＝（Q1－y）／IQR。如果L的絕對值超過給定的臨界值，則認為是異常值。臨界值的設定一般在1．5～3．0，通常認為取2．0附近的界值是合適的［5］。

穩健性回歸 假定某個胎齡x的新生兒出生體重y服從正態分布，穩健性回歸可用于識別數據中的異常值。常用的基于M估計的穩健性回歸通過最小化目標函數獲得未知參數的估計，其中n代表樣本量，ei表示殘差；特別地，ρ（e）＝e2對應參數的最小二乘估計。統計上，通過迭代再加權最小二乘方法獲得參數的估計值，每次迭代時殘差超過標準差的給定倍數的值將被識別為異常值，并從下次計算中剔除，直至估計的參數收斂。世界衛生組織（WHO）建議年齡別體重資料異常值的排除標準為4個標準差［8］。

高斯混合模型 針對數據雙峰的特點，假定某個胎齡x的新生兒出生體重y為兩個正態分布的復合分布f：一個為主要分布fp，另一個為次要分布fs，即f＝q·fp＋（1－q）·fs，其中，q代表樣本屬于這兩個正態分布的權重比例。這樣，觀測（x，y）屬于主要分布fp的概率為π＝q·fp／f。若π＞0．5，則認為（x，y）來自于主要分布，否則認為（x，y）為異常數據［3，7］。

3．曲線構建

獲得“干凈”數據后，直接使用各胎齡別數據的百分位數構建參考曲線比較粗劣。因此，常采用非參數方法或參數方法進行曲線平滑。非參數方法指各種統計平滑法，如三次樣條平滑［6］。最常用的參數方法莫過于Box－Cox指數－中位數－變異系數（lambdamedian－sigma，LMS）法［2，9，10］，并進一步推廣為基于位置參數、尺度參數和形狀參數的廣義可加模型（generalized additive model for location，scale and shape，GAMLSS）［6，11－12］。

（1）三次樣條方法 假定有n對測量值（xi，yi）（i＝1，…，n）且α≤x1≤x2≤…≤xn≤b，三次樣條函數f（x）定義為［a，b］上的函數，滿足：①f（x）在每個［xi，xi＋1］具有三次多項式形式f（x）＝β0i＋β1i（x－xi）＋β2i（x－xi）2＋β3i（x－xi）3；②f（x）及其前二階導數在xi連續。最小化懲罰性殘差平方和得到f（x）的懲罰性最小二乘估計，這里，λ（λ＞0）為給定的常數，在人體測量值擬合時取0．2左右是合適的［1］。

（2）LMS方法 LMS方法采用Box－Cox變換將胎齡別出生體重轉換為正態數值［10］。這個過程包括：估計每個胎齡t的三組模型參數，即正態轉換指數、體重中位數和變異系數；然后用三次樣條擬合三組參數獲得光滑的參考曲線。Cole和Green（1992）用Fisher得分法最大化懲罰性似然函數獲得光滑參數和模型參數的迭代估計。百分位數曲線表示為

Zα是標準正態分布的下百分位數。

（3）GAMLSS方法 GAMLSS方法允許對各種峰態非對稱分布進行建模并估計光滑的百分位數［11］。基于Box-Cox冪指數（Box-Cox power exponential，BCPE）分布的GAMLSS方法包含四個參數μ、σ、v和τ，分別對應胎齡別出生體重y的位置參數、尺度參數、偏度和峰度。GAMLSS對這些參數進行建模，g1（μ）＝h1（x），g2（σ）＝h2（x），g3（v）＝h3（x），g4（τ）＝h4（x）。其中，gi是對應參數的單調連接函數，如μ的恒等函數或對數函數；hi是關于協變量的參數或非參數函數，如時間t的線性函數或三次樣條函數。Rigby和Stasinopoulos（2005）采用Fisher得分算法最大化懲罰性似然函數進行光滑參數和模型參數的迭代估計，基于模型參數的估計值，可以用y的Z得分值獲得其任意百分位數估計。最近的研究建議使用對數連接函數對胎齡別出生體重構建百分位數曲線［12］。

4．擬合效果評價

采用赤池信息準則（Akaike information criteria，AIC）和殘差均方比較百分位數的擬合效果。其中，殘差均方定義為擬合的百分位數（P3，P10，P25，P50，P75，P90，P97）與實際值之差的平方和除以比較的百分位數的個數［2，15］。

5．統計軟件

本研究所有統計計算和建模過程均通過R軟件（Version 2．15．0，http：／／www．r－project．org）實現。

結 果

剔除缺失數據后，共納入167 288個初產單胎男孩的出生數據進行分析。從圖1看出，26～34周出生的新生兒，胎齡別出生體重呈明顯的雙峰分布或偏態分布。隨著胎齡的增加，雙峰分布逐漸消失，分布的對稱性有所改善。表1用統計數字特征描述了各個胎齡出生體重的分布情況。

圖1 單胎初產男嬰不同胎齡出生體重分布圖

以剔除盡可能少的異常點但又能捕獲分布的主要特征為基本依據，選取Tukey方法L統計量界值為2．0，高斯混合模型主要分布概率界值π為0．5，穩健性回歸界值為雙邊概率1／10000的3．89個標準差。最終，高斯混合模型在26～34周剔除18．9%（734／3885）的異常點，穩健性回歸在35～43周剔除0．9%（1493／163403）的異常點；而Tukey方法在26～34和35～43分別剔除5．9%（915／3885）和0．7%（1221／163403）的異常點。

圖2顯示了異常值剔除后胎齡別出生體重分布。其中，高斯混合模型有效地捕獲了26～34周出生的新生兒體重分布，Tukey方法在37～43周出生新生兒體重極端值的識別與穩健性回歸相近，均接近于剔除前的初始分位數。而從表1中35～36周出生新生兒體重的偏度和峰度系數來看，穩健性回歸更有效地剔除了由于孕周錯分導致的異常值。

圖3顯示了異常值剔除前后三次樣條光滑、LMS和GAMLSS三種方法擬合的胎齡別出生體重百分位數曲線。異常值剔除前，三種方法在35～43周的擬合結果相近，擬合值均接近于初始百分位數；但在26～34周，初始百分位數顯示小孕周出生體重較大孕周分散，提示存在系統誤差。異常值剔除后，三者擬合的曲線形狀相近，LMS和GAMLSS曲線的平滑效果優于三次樣條方法，尤其表現在26～27周樣本量較少的一段。表2比較了三種方法的擬合優度及用三種曲線劃分SGA和LGA的效果。AIC統計量及殘差均方均顯示GAMLSS方法的擬合效果與LMS方法相近，且優于三次樣條方法；但GAMLSS方法對SGA和LGA的判定更接近10%，因此擬合效果最好。比較圖3異常值剔除前后的擬合結果可見，異常值剔除主要糾正了26～34周各胎齡出生體重右尾部的分位數。

表1 三種異常值剔除方法胎齡別出生體重比較

圖2 三種異常數據剔除方法對胎齡別出生體重百分位數曲線構建的影響

表2 不同方法構建的參考曲線的擬合優度及SGA、AGA和LGA的分布

圖3 不同方法構建的參考曲線效果比較

討 論

在臨床應用中，出生體重參考曲線判別異常情況的能力備受重視。胎齡的錯分和出生體重極端值不僅影響參考曲線初始百分位數的估計，也會影響隨后的參考曲線平滑或建模［4，13］。本研究顯示，高斯混合模型較好地識別了多峰分布中的主要分布；Tukey方法和穩健性回歸對正態分布數據極端值的識別能力和效果比較接近；當數據呈現偏態和峰態分布時，穩健性回歸對異常數據的識別能力較好。因此，基于數據的分布特點和異常值的產生原因選用不同的方法，有利于改善異常數據的有效識別。

參考曲線構建時，利用三次樣條方法直接對分位數的經驗估計值進行平滑，方法直觀且對數據的分布假設較少，擬合結果接近經驗分位數，但效果并不理想，尤其是經驗百分位數曲線不夠光滑時，平滑效果較差［13］。采用懲罰性似然估計的LMS方法構建參考曲線已經得到廣泛驗證，即便出現隨時間變化的偏度系數，LMS方法得到的百分位數曲線也是可靠的。而且，LMS方法得到的L、M和S曲線能夠很好地反映觀測值的分布特征，便于更好地認識數據的深層結構。但是，如果觀測數據既呈偏態又呈峰態分布時，Box－Cox正態變換效果欠佳［11］。GAMLSS方法允許對中位數、標準差、偏度系數和峰度系數等四階矩進行建模，能夠極大地提高估計的精度；它的另一個優勢在于可以對胎齡、性別、產次等其他協變量進行建模，從而能在一個模型中使用全部數據，獲得樣本量較小分組的分布曲線的穩健估計［12］。

同其他統計建模一樣，生長曲線構建過程中參數的選擇至關重要。異常數據識別時，雖然沒有金標準，但剔除盡可能少的樣本而又能捕獲數據的基本特征成為廣泛采用的標準［3，6］；曲線構建的參數選擇，采用擬合優度統計量可方便地進行方法內和方法間擬合效果比較［9］。本文在異常值識別充分參考了文獻相關方法的參數選擇范圍，應用于本研究并進行異常值識別比例和剔除前后數據分布的比較。曲線構建時，我們也嘗試了0．1～0．4的其他光滑參數，平衡曲線光滑效果和擬合優度后選取了0．2；GAMLSS擬合時，我們也嘗試了恒等函數作為連接函數，但AIC統計量明顯比對數函數時大得多。最后，我們還比較了三種方法對SGA和LGA的識別能力，發現GAMLSS方法能夠識別9．58%的出生嬰兒為SGA，10．09%的為LGA，與國際上其他研究相近［14］。

新生兒出生體重僅僅是兒童生長發育指標的一個代表，本文闡述的參考曲線構建中常用的數據剔除和建模方法，可以推廣到其他生長曲線的構建，建模過程中統計方法的選擇和評價也可為其他相關工作提供參考。

［1］尚磊，徐勇勇，侯茹蘭，等．采用三次樣條函數擬合體重百分位數曲線．中國衛生統計，2001，18（5）：266-268．

［2］江梅．LMS法：一種適用建立肺功能全年齡段正常參考值曲線方法．中國衛生統計，2013，30（5）：766-768．

［3］He J，Xia H，Liu Y，et al．A new birthweight reference in Guangzhou，southern China，and its comparison with the global reference．Archives of diseases in childhood．2014：doi：10．1136／archdischild-2013-305923．

［4］Foster PJ，Kecojevic＇T．Reference Grow th Charts for Saudi Arabian Children and Adolescents．In：Vinod HD．Advances in Social Science Research Using R．New York：Springer，2010：107-128．

［5］Arbuckle TE，W ilkins R，Sherman GJ．Birth weight percentiles by gestational age in Canada．Obstetrics and gynecology，1993，81（1）：39-48．

［6］Bonellie S，Chalmers J，Gray R，et al．Centile charts for birthweight for gestational age for Scottish singleton births．BMC pregnancy and childbirth，2008，8：5．

［7］Liu Z，Zhang J，Zhao B，et al．Population-based reference for birth weight for gestational age in northern China．Early human development，2014，90（4）：177-187．

［8］Centers for Disease Controland Prevention．Cut-offs to define outliers in the 2000 CDC Grow th Charts．Available from：http：／／www．cdc．gov／nccdphp／dnpa／grow thcharts／resources／BIV-cutoffs．pdf

［9］Flegal KM．Curve smoothing and transformations in the development of grow th curves．The American journal of clinical nutrition，1999，（1）：163S-165S．

［10］Cole TJ，Green PJ．Smoothing reference centile curves：the LMS method and penalized likelihood．Statistics in medicine．1992，11（10）：1305-1319．

［11］Rigby RA，Stasinopoulos DM．Generalized additivemodels for location，scale and shape．Applied statistics，2005，54：507-554．

［12］Cole TJ，Statnikov Y，Santhakumaran S，et al．Birth weight and longitudinal grow th in infants born below 32 weeks’gestation：a UK population study．Archives of disease in childhood-fetal and neonatal edition，2014，99：34-40．

［13］Ramos F，Perez G，Jane M，et al．Construction of the birth weight by gestational age population reference curves of Catalonia（Spain）：Methods and development．Gaceta santitaria，2009，23：76-81．

［14］Olsen IE，Groveman SA，Lawson ML，et al．New intrauterine grow th curves based on United States data．Pediatrics，2010，125：e214-224

（責任編輯：郭海強）

Com parison of Statistical Methods for the Construction of Grow th Reference Curves

Liu Yu，He Jianrong，Guo Yong，et al．（GuangzhouWomenandChildren’sMedicalCenter（510623），Guangzhou）

ObjectiveTo find out amethod to construct grow th curves for various anthropology measurements．MethodsIllustrated by the construction of a new birth weight reference for gestational age in Guangzhou，the best data pre-processing procedure was determined for outlier detection and removal after comparing with 3 differentmethods：Tukey’smethod，robust regression and Gaussian mixture model．Then cubic spline，LMS and GAMLSS were compared for the construction of grow th curves．ResultsGaussianm ixturemodel had an outstanding performance in the detection of outliers for datawithmultimodal distribution and compared with Tukey’smethod，robust regression seemedmore reliable for the unimodal asymmetric data．In the construction of centile charts，GAMLSSmodel gavemostaccurate estimate compared with cubic spline and LMSmethod．ConclusionMethods based on the distribution characteristics should be used for data pre-processing．The fourthmoment of the distribution can bemodeled by GAMLSS，so grow th curves fitted by GAMLSS aremuch smoother and the fitting errors are much fewer．

Grow th curve；Tukey’smethod；Robust regression；Gaussian mixturemodel；LMS；GAMLSS

*：廣州市科技計劃重大科技專項（2011Y2－00025）；廣州市科技計劃國際科技交流與合作專項（2012J5100038）

1．廣州市婦女兒童醫療中心（510623）

2．中山大學數學與計算科學學院統計系

△通訊作者：邱琇，E-mail：qxiu0161＠163．com