簡(jiǎn)介面向2022年全國(guó)大學(xué)生數(shù)學(xué)建模競(jìng)賽C題備賽者這份PDF完整收錄了古代玻璃文物成分分析與鑒別問(wèn)題的解題方案與配套代碼說(shuō)明聚焦數(shù)學(xué)建模在考古數(shù)據(jù)場(chǎng)景下的實(shí)際應(yīng)用。內(nèi)容圍繞數(shù)據(jù)預(yù)處理、主成分分析、因子分析、敏感性分析及灰色關(guān)聯(lián)度分析等核心方法展開(kāi)系統(tǒng)覆蓋了從數(shù)據(jù)探索、風(fēng)化系數(shù)定義、亞類劃分、分類預(yù)測(cè)到不同玻璃類型化學(xué)成分關(guān)聯(lián)性比較的全流程建模思路并附有關(guān)鍵代碼與結(jié)果解釋便于讀者對(duì)照理解。資源共1個(gè)文件大小3.55MBPDF格式方便在線閱讀與打印。目前已有927人學(xué)習(xí)下載其問(wèn)題拆解和模型選擇邏輯具有較高參考價(jià)值尤其適合備賽數(shù)學(xué)建模、需要參考完整賽題分析框架的高校學(xué)生與建模愛(ài)好者。讀者可快速掌握風(fēng)化數(shù)據(jù)場(chǎng)景下的降維、因子轉(zhuǎn)換、魯棒性檢驗(yàn)與灰色關(guān)聯(lián)度建模技巧并將其遷移到其他成分分析類問(wèn)題中。1. 風(fēng)化系數(shù)為什么風(fēng)化效應(yīng)會(huì)讓成分?jǐn)?shù)據(jù)說(shuō)反話古代玻璃制品的成分分析最反直覺(jué)的地方在于風(fēng)化并不是簡(jiǎn)單地讓某些成分流失而是讓整組數(shù)據(jù)都在說(shuō)謊。同樣是埋藏環(huán)境下高鉀玻璃風(fēng)化后二氧化硅占比顯著升高、氧化鉀顯著降低鉛鋇玻璃卻完全反過(guò)來(lái)二氧化硅占比顯著降低、氧化鉛和氧化鋇顯著升高。如果拿著風(fēng)化后的檢測(cè)數(shù)據(jù)直接做分類或回歸模型學(xué)到的很可能是風(fēng)化程度而不是玻璃類型本身。這份 2022 年數(shù)學(xué)建模 C 題的標(biāo)準(zhǔn)解法核心思路就是先定義風(fēng)化系數(shù)把風(fēng)化程度量化再用主成分分析降維、因子分析反推未風(fēng)化成分、灰色關(guān)聯(lián)度分析化學(xué)成分關(guān)聯(lián)——這一整套鏈路解決的不只是考古問(wèn)題更是觀測(cè)值被環(huán)境篡改后如何建模的通用范式。2. 數(shù)據(jù)探索先行三個(gè)反直覺(jué)事實(shí)決定建模方向建模比賽里最常見(jiàn)的翻車點(diǎn)不是模型不夠高級(jí)而是沒(méi)看數(shù)據(jù)就套模型。這道題最值得借鑒的地方恰恰是把數(shù)據(jù)探索放在了建模之前而且探索結(jié)論直接決定了后續(xù)每個(gè)模型的選擇。2.1 兩類玻璃的風(fēng)化方向截然相反玻璃的主要成分是二氧化硅但為了降低熔融溫度古人會(huì)添加助熔劑。加入鉛礦石、鋇礦石得到的是鉛鋇玻璃加入草木灰得到的是高鉀玻璃。這些助熔劑殘留物成為玻璃的次要成分氧化鉛PbO、氧化鋇BaO、氧化鉀K2O、氧化鈣CaO、氧化鎂MgO、氧化鋁Al2O3等。風(fēng)化效應(yīng)是玻璃與埋藏環(huán)境長(zhǎng)期進(jìn)行化學(xué)元素交換的結(jié)果。關(guān)鍵在于兩類玻璃的風(fēng)化產(chǎn)物方向完全相反成分變化方向鉛鋇玻璃高鉀玻璃二氧化硅SiO2顯著降低顯著升高氧化鉛 / 氧化鋇顯著升高變化不顯著氧化鉀變化不顯著顯著降低這個(gè)結(jié)論來(lái)自文獻(xiàn)中對(duì)鉛鋇玻璃、鉀玻璃風(fēng)化機(jī)理的說(shuō)明也在附件表單中得到了驗(yàn)證。高鉀玻璃風(fēng)化后表面富硅層形成二氧化硅相對(duì)占比被抬高鉛鋇玻璃風(fēng)化后鉛、鋇成分相對(duì)富集二氧化硅被稀釋。同一套風(fēng)化機(jī)制作用在兩類玻璃上產(chǎn)生了方向相反的表觀變化。這帶來(lái)一個(gè)直接后果如果拿風(fēng)化后的原始成分?jǐn)?shù)據(jù)做聚類高鉀風(fēng)化樣本和鉛鋇未風(fēng)化樣本可能在某些維度上靠得很近因?yàn)楦哜涳L(fēng)化后 SiO2 升高、K2O 降低鉛鋇未風(fēng)化時(shí)恰好也是 SiO2 高、鉛鋇低。數(shù)據(jù)探索的價(jià)值在這里就體現(xiàn)出來(lái)了——它告訴你必須先給風(fēng)化程度一個(gè)定量指標(biāo)否則后續(xù)所有分類都建立在流沙上。2.2 風(fēng)化系數(shù)的定義邏輯與閾值判定論文里給出了兩個(gè)風(fēng)化系數(shù)定義這是整個(gè)建模的地基。對(duì)鉛鋇玻璃風(fēng)化系數(shù)定義為# 鉛鋇玻璃風(fēng)化導(dǎo)致 PbO、BaO 升高SiO2 降低 # 所以比值越大風(fēng)化越嚴(yán)重 lead_barium_weathering (PbO BaO) / SiO2對(duì)高鉀玻璃風(fēng)化系數(shù)定義為# 高鉀玻璃風(fēng)化導(dǎo)致 SiO2 升高K2O、CaO、MgO、Al2O3、Fe2O3、CuO 降低 # 所以 SiO2 在分子上其他成分在分母上 potassium_weathering SiO2 / (K2O CaO MgO Al2O3 Fe2O3 CuO)這兩個(gè)定義不是拍腦袋想出來(lái)的而是嚴(yán)格對(duì)應(yīng)風(fēng)化方向風(fēng)化系數(shù)應(yīng)該和風(fēng)化程度單調(diào)相關(guān)。鉛鋇玻璃風(fēng)化越重PbO、BaO 越多SiO2 越少比值自然越大高鉀玻璃風(fēng)化越重SiO2 越多其他成分越少比值也越大。有了風(fēng)化系數(shù)之后閾值判定就有了依據(jù)。觀察數(shù)據(jù)分布可以發(fā)現(xiàn)一個(gè)典型現(xiàn)象鉛鋇玻璃-A 紋類型的風(fēng)化系數(shù)在 0.879 和 1.3 之間存在明顯突變中間幾乎沒(méi)有過(guò)渡值高鉀玻璃的風(fēng)化系數(shù)在 7.23 和 12.41 之間也存在跳躍。這種斷層說(shuō)明風(fēng)化樣本和未風(fēng)化樣本在定量指標(biāo)上是可分的于是把閾值分別定在 1.0 和 10.0。這里有個(gè)容易被忽略的細(xì)節(jié)閾值不是從統(tǒng)計(jì)檢驗(yàn)里算出來(lái)的而是直接觀察數(shù)據(jù)突變得到的。這個(gè)做法在數(shù)學(xué)建模里是合理的——你是在找一個(gè)能把樣本分成兩堆的分界點(diǎn)而不是在估計(jì)一個(gè)總體參數(shù)。不過(guò)論文里也承認(rèn)存在少量標(biāo)簽不一致的情況比如某樣本的標(biāo)注與風(fēng)化系數(shù)判別的結(jié)果矛盾最后以風(fēng)化系數(shù)的判定為準(zhǔn)做了手工調(diào)整。這種做法在比賽里無(wú)可厚非但放到真實(shí)考古場(chǎng)景中就意味著你需要第三方證據(jù)來(lái)驗(yàn)證標(biāo)注錯(cuò)了而不是風(fēng)化系數(shù)定義錯(cuò)了。2.3 數(shù)據(jù)探索發(fā)現(xiàn)的其他規(guī)律除了風(fēng)化方向相反之外數(shù)據(jù)探索還發(fā)現(xiàn)了三個(gè)具體現(xiàn)象第一紋飾 B 的玻璃文物全部是高鉀玻璃、全部風(fēng)化、顏色全部是藍(lán)綠色。這個(gè)組合的高一致性說(shuō)明紋飾和產(chǎn)地、時(shí)代高度相關(guān)。第二紋飾 A 的高鉀玻璃全部無(wú)風(fēng)化而鉛鋇玻璃風(fēng)化與未風(fēng)化比例約為 11:5。第三紋飾 C 的高鉀玻璃全部無(wú)風(fēng)化鉛鋇玻璃風(fēng)化與未風(fēng)化比例為 17:73。把這三個(gè)現(xiàn)象放在一起可以得到兩個(gè)推論同等條件下鉛鋇玻璃比高鉀玻璃更容易風(fēng)化紋飾 B 所代表的環(huán)境對(duì)玻璃有更強(qiáng)的風(fēng)化效應(yīng)因?yàn)檫B相對(duì)耐風(fēng)化的高鉀玻璃都全部風(fēng)化了。這些結(jié)論直接在論文中作為問(wèn)題一的回答基礎(chǔ)也為問(wèn)題二的分類模型提供了先驗(yàn)知識(shí)分類時(shí)必須把是否風(fēng)化作為分層變量先區(qū)分風(fēng)化與未風(fēng)化再做類型判別。這就是為什么后面的主成分分析要分四類數(shù)據(jù)分別做——高鉀風(fēng)化、高鉀未風(fēng)化、鉛鋇風(fēng)化、鉛鋇未風(fēng)化——而不是把所有樣本丟進(jìn)一個(gè) PCA 里。如果混在一起降維主成分會(huì)被風(fēng)化差異主導(dǎo)類型差異反而被淹沒(méi)。3. 主成分與因子分析把風(fēng)化后的成分反推回未風(fēng)化題目最硬核的問(wèn)題一是預(yù)測(cè)風(fēng)化前的成分含量。這本質(zhì)上是一個(gè)反演問(wèn)題你只能測(cè)到風(fēng)化后的成分但想還原風(fēng)化前的成分。這個(gè)問(wèn)題的解決思路非常巧妙——用因子分析建模風(fēng)化過(guò)程再對(duì)負(fù)載矩陣求逆反推。3.1 為什么必須先做主成分分析降維先看數(shù)據(jù)維度。表單 2 中每個(gè)樣本有 14 個(gè)化學(xué)成分變量而樣本量只有幾十個(gè)。直接用 14 維數(shù)據(jù)做因子分析會(huì)遇到兩個(gè)問(wèn)題一是運(yùn)算量偏大二是很多成分之間高度相關(guān)存在多重共線性直接做回歸或因子分析結(jié)果會(huì)被干擾成分帶偏。主成分分析的作用就是把 14 個(gè)原始變量壓縮成少數(shù)幾個(gè)互不相關(guān)的主成分。每個(gè)主成分是原始變量的線性組合第一個(gè)主成分承載最大方差第二個(gè)次之依此類推。論文里的具體做法是按玻璃類型和風(fēng)化狀態(tài)分成四類數(shù)據(jù)對(duì)每一類分別做 PCA選取累計(jì)貢獻(xiàn)率達(dá)到 85% 以上的主成分。以鉛鋇未風(fēng)化類型為例前 6 個(gè)主成分的累計(jì)貢獻(xiàn)率達(dá)到 88.89%對(duì)應(yīng)的主要化學(xué)成分是 SiO2、Na2O、K2O、CaO、MgO、Al2O3。高鉀未風(fēng)化類型的前 5 個(gè)成分貢獻(xiàn)率同樣達(dá)到 88.89%。這意味著原始 14 維信息被壓縮到 5-8 維同時(shí)保留了絕大部分方差。PCA 的具體實(shí)現(xiàn)并不復(fù)雜import numpy as np def pca_components(X, n_componentsNone, threshold0.85): # X: 樣本矩陣每行一個(gè)樣本每列一個(gè)化學(xué)成分 # 第一步標(biāo)準(zhǔn)化 X_mean X.mean(axis0) X_std X.std(axis0) X_norm (X - X_mean) / X_std # 第二步計(jì)算相關(guān)系數(shù)矩陣 # 注意這里用的是相關(guān)系數(shù)矩陣而非協(xié)方差矩陣 # 因?yàn)椴煌瘜W(xué)成分量綱相同但方差差異大標(biāo)準(zhǔn)化后更公平 R np.corrcoef(X_norm, rowvarFalse) # 第三步特征值分解 eigenvalues, eigenvectors np.linalg.eigh(R) idx np.argsort(eigenvalues)[::-1] eigenvalues eigenvalues[idx] eigenvectors eigenvectors[:, idx] # 第四步按累計(jì)貢獻(xiàn)率選主成分個(gè)數(shù) total_var eigenvalues.sum() cum_ratio np.cumsum(eigenvalues) / total_var if n_components is None: n_components np.searchsorted(cum_ratio, threshold) 1 # 主成分得分標(biāo)準(zhǔn)化后的數(shù)據(jù)投影到特征向量上 scores X_norm eigenvectors[:, :n_components] return scores, eigenvectors[:, :n_components], cum_ratio這里有個(gè)關(guān)鍵細(xì)節(jié)做 PCA 之前必須先標(biāo)準(zhǔn)化。如果直接對(duì)原始成分?jǐn)?shù)據(jù)做特征值分解SiO2 這種含量經(jīng)常在 60% 以上的變量會(huì)主導(dǎo)第一主成分而 PbO 這種含量波動(dòng)大的變量卻因?yàn)閿?shù)值小被忽略。標(biāo)準(zhǔn)化之后每個(gè)變量在降維決策中的權(quán)重才真正由其相關(guān)性決定。3.2 因子分析模型的核心邏輯因子分析假設(shè)可觀測(cè)的化學(xué)成分向量 x 由公共因子 f 線性組合而成x μ A f ε其中 μ 是成分期望向量f 是公共因子向量這里對(duì)應(yīng)未風(fēng)化時(shí)的成分特征A 是負(fù)載矩陣ε 是特殊因子即隨機(jī)誤差。這個(gè)模型的物理含義是風(fēng)化的過(guò)程是連續(xù)的、由多個(gè)公共因子驅(qū)動(dòng)的某個(gè)時(shí)刻觀測(cè)到的成分比例是原始成分經(jīng)過(guò)風(fēng)化演化后在當(dāng)前時(shí)刻的表現(xiàn)。公共因子 f 數(shù)量小于原始變量數(shù)量這就是降維的本質(zhì)。關(guān)鍵的建模步驟在于如何獲取負(fù)載矩陣 A。論文中的做法是對(duì)標(biāo)準(zhǔn)化后數(shù)據(jù)的相關(guān)系數(shù)矩陣做特征值分解取前 m 個(gè)特征值 λ 和對(duì)應(yīng)的特征向量 u構(gòu)造A [sqrt(λ1)·u1, sqrt(λ2)·u2, ..., sqrt(λm)·um]也就是把特征值開(kāi)方后乘到特征向量上。這個(gè) A 的列就對(duì)應(yīng)各個(gè)公共因子在各個(gè)原始變量上的載荷。獲取公共因子的表達(dá)式需要對(duì)式 (3) 做變換。把 x μ A f ε 稍作整理得到 f 關(guān)于 x 的線性表示f A^{-1} (x - μ)這就是由風(fēng)化后成分反推未風(fēng)化成分的轉(zhuǎn)換公式。論文中給出的是鉛鋇玻璃-C 紋類型的轉(zhuǎn)換系數(shù)矩陣——一個(gè) 8 行 10 列的系數(shù)表每一行對(duì)應(yīng)一個(gè)公共因子每一列對(duì)應(yīng)一個(gè)原始化學(xué)成分。把風(fēng)化后的實(shí)測(cè)成分代入就能得到未風(fēng)化時(shí)的預(yù)測(cè)值。用代碼來(lái)演示這個(gè)反推過(guò)程會(huì)非常直觀import numpy as np def factor_reverse(observed, mu, A): 由風(fēng)化后成分反推未風(fēng)化成分 observed: 風(fēng)化后檢測(cè)到的化學(xué)成分向量 mu: 訓(xùn)練集成分均值向量 A: 負(fù)載矩陣shape (n_components, n_features) # 中心化減去均值 centered observed - mu # 求負(fù)載矩陣的偽逆 # 注意A 不一定是方陣主成分?jǐn)?shù) 特征數(shù) # 所以不能直接求逆要用 pinv A_inv np.linalg.pinv(A) # 反推公共因子值即未風(fēng)化時(shí)的主成分得分 factors A_inv centered # 再由公共因子重構(gòu)原始空間 recovered mu A.T factors return recovered這段代碼里最值得注意的是np.linalg.pinv。比賽復(fù)盤時(shí)很多人會(huì)在這里踩坑負(fù)載矩陣 A 的形狀是 (主成分?jǐn)?shù), 特征數(shù))不是方陣直接調(diào)用np.linalg.inv會(huì)直接報(bào)錯(cuò)必須用偽逆。偽逆的本質(zhì)是最小二乘意義下的最優(yōu)解當(dāng) A 的列數(shù)多于行數(shù)時(shí)它給出的解是范數(shù)最小的那一個(gè)。3.3 主成分個(gè)數(shù)與累計(jì)貢獻(xiàn)率的選擇論文中針對(duì)鉛鋇玻璃-C 紋選取了 8 個(gè)主成分保留 97% 以上的成分構(gòu)成信息。這個(gè)選擇不是越高越好而是有講究的。主成分個(gè)數(shù)太少會(huì)丟失重要的成分變化信息太多則模型復(fù)雜度和過(guò)擬合風(fēng)險(xiǎn)上升。數(shù)學(xué)建模中通常以 85% 作為累計(jì)貢獻(xiàn)率的下限但對(duì)這道題來(lái)說(shuō)由于最終要做因子反演主成分個(gè)數(shù)直接影響負(fù)載矩陣 A 的列數(shù)進(jìn)而影響反推結(jié)果的穩(wěn)定性。這里涉及一個(gè)數(shù)字上的權(quán)衡如果只選前 3 個(gè)主成分累計(jì)貢獻(xiàn)率可能只有 70%反推出來(lái)的未風(fēng)化成分會(huì)嚴(yán)重失真如果選到 8 個(gè)貢獻(xiàn)率到了 97%但 A 的維度變大偽逆計(jì)算對(duì)噪聲的放大效應(yīng)也會(huì)增強(qiáng)。論文里的選擇是偏保守的優(yōu)先保證還原精度犧牲了一部分運(yùn)算效率。從論證角度看這個(gè)選擇是站得住腳的。因子分析在這里不是一個(gè)預(yù)測(cè)精確數(shù)值的工具而是一個(gè)還原成分比例的工具。對(duì)于這種目標(biāo)盡可能保留信息的優(yōu)先級(jí)應(yīng)該高于模型簡(jiǎn)潔性。這也解釋了為什么問(wèn)題三的分類只用區(qū)間判斷就夠了——分類關(guān)心的是樣本落在哪個(gè)范圍不需要精確還原而問(wèn)題一需要精確的成分比例必須用因子分析。4. 分類、亞類劃分與敏感性分析建模落地與常見(jiàn)問(wèn)題排查問(wèn)題二、三、四本質(zhì)上是一套分類 亞類劃分 穩(wěn)定性檢驗(yàn)的組合拳。這部分看起來(lái)沒(méi)有因子分析那么驚艷但恰恰是建模完整性的體現(xiàn)也是數(shù)據(jù)挖掘類賽題最容易丟分的地方。4.1 亞類劃分為什么用最短距離法聚類二氧化硅問(wèn)題二要求在每一類高鉀/鉛鋇、風(fēng)化/未風(fēng)化內(nèi)部再做亞類劃分。論文的做法分兩步第一步用主成分分析找出貢獻(xiàn)率最高的成分。四類數(shù)據(jù)各自的 PCA 結(jié)果中SiO2 都是差異最大的一個(gè)成分。于是選擇二氧化硅作為亞類劃分的依據(jù)。第二步用最短距離法single linkage對(duì)二氧化硅數(shù)據(jù)聚類。最短距離法的核心定義是兩類之間的距離等于兩類中最近的兩個(gè)樣本之間的距離。D_pq min { d(x_i, x_j) : x_i ∈ G_p, x_j ∈ G_q }合并后新類到其他類的距離遞推公式為D_rk min { D_pk, D_qk }這個(gè)方法的優(yōu)點(diǎn)是計(jì)算簡(jiǎn)單每次只要找距離矩陣中的最小非對(duì)角元素即可合并。缺點(diǎn)在算法層面也很明顯它對(duì)噪聲和離群點(diǎn)非常敏感只要有一個(gè)異常樣本離得遠(yuǎn)就可能讓兩個(gè)本不該合并的類被拉在一起。論文實(shí)際是把二氧化硅數(shù)據(jù)聚類成兩個(gè)聚點(diǎn)然后取兩個(gè)聚點(diǎn)的平均值作為分界把數(shù)據(jù)分成兩個(gè)區(qū)間。以鉛鋇未風(fēng)化類型為例兩個(gè)聚點(diǎn)分別是 64.82 和 52.49分界點(diǎn)在兩者平均值附近最終得到區(qū)間 [50.61, 58.70) 和 [58.70, 69.71]。這種區(qū)間劃分方法在數(shù)學(xué)建模中的合理性在于它的結(jié)果可解釋性極強(qiáng)。給定一個(gè)新樣本只需要檢查它的 SiO2 落在哪個(gè)區(qū)間就能直接判定亞類。不需要重新訓(xùn)練模型不需要調(diào)參一個(gè) Excel 公式就能完成分類。4.2 敏感性分析的雙面結(jié)論說(shuō)明了什么敏感性分析是這道題的一個(gè)關(guān)鍵得分點(diǎn)論文中給出了兩個(gè)看似矛盾的結(jié)論問(wèn)題二對(duì)四類已知數(shù)據(jù)的二氧化硅添加 3% 隨機(jī)擾動(dòng)后再分類發(fā)現(xiàn)有三類數(shù)據(jù)的分類誤差超過(guò) 20%于是判定分類結(jié)果的敏感性較強(qiáng)。但問(wèn)題三對(duì)表單 3 的未知數(shù)據(jù)做同樣的 3% 擾動(dòng)分類結(jié)果與原來(lái)完全一致判定敏感性不強(qiáng)。同一個(gè)方法出現(xiàn)兩種結(jié)論這不是模型出錯(cuò)而是兩類數(shù)據(jù)本身的分布特性不同。訓(xùn)練集數(shù)據(jù)表單 2的二氧化硅值分布比較分散很多樣本靠近區(qū)間邊界3% 的擾動(dòng)足以讓它們?cè)竭^(guò)邊界跳到另一個(gè)區(qū)間而表單 3 的未知樣本則大多落在區(qū)間內(nèi)部距離邊界足夠遠(yuǎn)同樣幅度的擾動(dòng)不會(huì)改變歸屬。這里有一個(gè)值得深挖的細(xì)節(jié)區(qū)間分類法天然帶有邊界敏感問(wèn)題。只要樣本落在邊界附近微小擾動(dòng)就會(huì)導(dǎo)致分類翻轉(zhuǎn)。解決方法不是改用更復(fù)雜的模型而是檢查樣本在區(qū)間內(nèi)的位置。如果所有樣本都離邊界較遠(yuǎn)說(shuō)明分類結(jié)果是穩(wěn)定的如果存在大量靠近邊界的樣本則需要考慮把概率輸出而不是硬分類作為結(jié)果。4.3 灰色關(guān)聯(lián)度分析的完整計(jì)算流程問(wèn)題四要求分析化學(xué)成分之間的關(guān)聯(lián)關(guān)系論文選用灰色關(guān)聯(lián)度分析。這個(gè)方法特別適合考古數(shù)據(jù)因?yàn)樗灰髷?shù)據(jù)滿足正態(tài)分布或線性關(guān)系假設(shè)在樣本量小、數(shù)據(jù)不完全的情況下也能給出穩(wěn)定的相似性度量。計(jì)算流程分為五步import numpy as np def grey_relational_analysis(reference, comparison, rho0.5): reference: 參考數(shù)列形狀 (n,) comparison: 比較數(shù)列形狀 (m, n)每行是一個(gè)化學(xué)成分序列 rho: 分辨系數(shù)一般取 0.5 # 第一步無(wú)量綱化均值化 ref reference / reference.mean() comp comparison / comparison.mean(axis1, keepdimsTrue) # 第二步計(jì)算差值序列 delta np.abs(comp - ref) # 形狀 (m, n) # 第三步求兩級(jí)最小差和最大差 min_delta delta.min() max_delta delta.max() # 第四步計(jì)算關(guān)聯(lián)系數(shù) # 公式: ksi (min rho*max) / (delta rho*max) ksi (min_delta rho * max_delta) / (delta rho * max_delta) # 第五步關(guān)聯(lián)度 關(guān)聯(lián)系數(shù)均值 degree ksi.mean(axis1) return degree關(guān)于分辨系數(shù) ρ論文特意提到 ρ ≤ 0.5463 時(shí)分辨力最好一般取 0.5。ρ 越小關(guān)聯(lián)系數(shù)之間的差異越明顯分辨力越強(qiáng)ρ 越大所有系數(shù)向 1 靠攏區(qū)分度降低。取 0.5 是權(quán)衡分辨力和穩(wěn)定性的常見(jiàn)選擇。對(duì)四類玻璃分別計(jì)算關(guān)聯(lián)度并排序論文得到的結(jié)果是玻璃類型化學(xué)成分關(guān)聯(lián)度排序鉛鋇未風(fēng)化SiO2 CaO K2O Al2O3 MgO Na2O高鉀未風(fēng)化SiO2 CaO K2O MgO Na2O鉛鋇風(fēng)化SiO2 Na2O K2O CaO MgO高鉀風(fēng)化SiO2 Na2O K2O四類數(shù)據(jù)的共同點(diǎn)是 SiO2 都是關(guān)聯(lián)度最高的成分這符合玻璃以二氧化硅為主體結(jié)構(gòu)的基本事實(shí)。差異在于次要成分的關(guān)聯(lián)順序不同未風(fēng)化樣本中 CaO、K2O 與主體成分的關(guān)聯(lián)更強(qiáng)而風(fēng)化樣本中 Na2O 的關(guān)聯(lián)排名上升。這說(shuō)明風(fēng)化過(guò)程確實(shí)改變了化學(xué)成分之間的耦合關(guān)系也側(cè)面印證了風(fēng)化系數(shù)定義的合理性。4.4 避坑指南數(shù)學(xué)建模中最容易翻車的五個(gè)細(xì)節(jié)拿這套題復(fù)盤以下五個(gè)坑幾乎每年都有參賽隊(duì)踩進(jìn)去?？右恢苯訉?duì)全部樣本做聚類不分風(fēng)化狀態(tài)?，F(xiàn)象聚類結(jié)果總是把風(fēng)化樣本單獨(dú)分出來(lái)類型反而混在一起。原因風(fēng)化差異的方差遠(yuǎn)大于類型差異無(wú)監(jiān)督算法優(yōu)先捕捉方差最大的方向。解決先按類型 是否風(fēng)化分成四個(gè)子集再分別建模。這套題里所有分析都是四類分開(kāi)做的這是最大的結(jié)構(gòu)性決策。坑二風(fēng)化系數(shù)定義方向錯(cuò)了?，F(xiàn)象系數(shù)大的樣本反而顯示未風(fēng)化。原因把分子分母放反了。鉛鋇玻璃風(fēng)化后是 PbO、BaO 增加所以要把它們放分子高鉀玻璃風(fēng)化后是 SiO2 增加所以要把 SiO2 放分子。解決定義系數(shù)后先拿已知標(biāo)簽的樣本做驗(yàn)證如果系數(shù)和標(biāo)簽明顯矛盾優(yōu)先查定義而不是查數(shù)據(jù)?？尤蜃臃治龇赐茣r(shí)用inv而不是pinv。現(xiàn)象程序報(bào)錯(cuò)說(shuō)矩陣不是方陣。原因負(fù)載矩陣行數(shù)是主成分?jǐn)?shù)列數(shù)是特征數(shù)兩者不相等。解決改用np.linalg.pinv并在注釋里寫清楚是最小二乘意義下的逆。坑四敏感性分析只看對(duì)不對(duì)不看差多遠(yuǎn)?，F(xiàn)象隨機(jī)擾動(dòng)后分類結(jié)果變了就斷言模型不穩(wěn)。原因沒(méi)有區(qū)分樣本是落在區(qū)間中心還是區(qū)間邊界。解決除了報(bào)告分類準(zhǔn)確率還要報(bào)告每個(gè)樣本到最近區(qū)間邊界的距離分布。這是驗(yàn)證區(qū)間分類法可靠性的最直接指標(biāo)?？游寤疑P(guān)聯(lián)度分析沒(méi)有做無(wú)量綱化?，F(xiàn)象關(guān)聯(lián)度排序完全由數(shù)值大的成分主導(dǎo)SiO2 永遠(yuǎn)是第一。原因SiO2 的含量是其他成分的幾十倍絕對(duì)差被它壟斷。解決必須先做初值化或均值化處理把每個(gè)序列拉到同一量級(jí)否則整個(gè)分析沒(méi)有意義。5. 設(shè)計(jì)取舍這套方法為什么能work又會(huì)在哪里失效把整套建模鏈路拆開(kāi)看最值得思考的設(shè)計(jì)取舍是先定義風(fēng)化系數(shù)再走因子分析這兩步的配合。風(fēng)化系數(shù)本質(zhì)上是一個(gè)領(lǐng)域知識(shí)驅(qū)動(dòng)的特征工程它把物理風(fēng)化機(jī)理壓縮成一個(gè)標(biāo)量。它的作用是制造一個(gè)分水嶺讓后續(xù)所有模型在高信噪比的子空間里工作。這是整道題最關(guān)鍵的一步因?yàn)轱L(fēng)化方向相反的數(shù)據(jù)混在一起時(shí)任何純數(shù)據(jù)驅(qū)動(dòng)的算法都會(huì)被誤導(dǎo)。因子分析的設(shè)計(jì)則是把風(fēng)化建模成線性變換然后靠矩陣求逆反推。這里面隱含了一個(gè)強(qiáng)假設(shè)風(fēng)化過(guò)程在嚴(yán)重風(fēng)化之前是近似線性的。論文的問(wèn)題假設(shè)部分明確寫了這一點(diǎn)。這個(gè)假設(shè)在考古場(chǎng)景中是否成立短期內(nèi)規(guī)律是成立的但長(zhǎng)期風(fēng)化會(huì)有非線性拐點(diǎn)比如表面富硅層形成后反應(yīng)速率會(huì)明顯下降。所以這套方法對(duì)中度風(fēng)化樣本的還原效果好對(duì)嚴(yán)重風(fēng)化的樣本會(huì)系統(tǒng)性偏差。再看區(qū)間分類法。它的優(yōu)缺點(diǎn)在敏感性分析中暴露得非常充分好處是透明任何人都能看懂分類規(guī)則壞處是邊界處沒(méi)有概率概念樣本落在 [58.69, 58.71] 和 [60.00, 60.01] 被同等對(duì)待后者明明更穩(wěn)定。如果換成高斯混合模型或核密度估計(jì)可以給出更平滑的分類概率但解釋成本會(huì)上升。在數(shù)學(xué)建模競(jìng)賽的評(píng)分體系里可解釋性往往比精確度更值錢這個(gè)取舍是正確的?；疑P(guān)聯(lián)度分析是整套方法里最輕的一環(huán)它不需要分布假設(shè)也不存在過(guò)擬合問(wèn)題計(jì)算量極小。它的局限在于關(guān)聯(lián)度是對(duì)幾何形狀相似性的度量不一定是因果關(guān)系的度量。論文把關(guān)聯(lián)度排序差異作為風(fēng)化影響成分關(guān)聯(lián)的證據(jù)這個(gè)推斷方向只能說(shuō)合理不能說(shuō)嚴(yán)格。從整體評(píng)價(jià)這套方法它贏在流程的閉環(huán)。數(shù)據(jù)探索發(fā)現(xiàn)風(fēng)化方向定義系數(shù)量化風(fēng)化程度PCA 降維去噪因子分析反演區(qū)間分類加敏感性檢驗(yàn)灰色關(guān)聯(lián)度做關(guān)系分析每一步的輸入輸出都銜接得嚴(yán)絲合縫。它不依賴任何黑匣子模型每一步都能被檢驗(yàn)、被解釋、被復(fù)現(xiàn)。正是這種透明 可驗(yàn)證的特質(zhì)讓它成為數(shù)學(xué)建模題的一個(gè)典范解。我自己的一個(gè)習(xí)慣是拿到任何一份成分?jǐn)?shù)據(jù)或觀測(cè)數(shù)據(jù)第一件事不是跑 PCA而是先定義好什么是一個(gè)樣本的本質(zhì)身份什么又是被環(huán)境污染的噪聲。從那以后我每次做這類反演任務(wù)都強(qiáng)制先走一遍定義污染系數(shù) → 驗(yàn)證方向 → 再看模型這個(gè)流程省下的返工時(shí)間遠(yuǎn)超想象。希望幫到你。本文還有配套的精品資源點(diǎn)擊獲取