馬 由,湯 艷,解 斐
(中國電子科技集團第十五研究所 軟件測評中心,北京 100083)
本文采用兩種模型來實現(xiàn)軟件缺陷數(shù)預測,即一般的線性回歸模型和貝葉斯網絡模型,并基于LASSO嵌入式特征選擇算法選擇最優(yōu)的因子集合。LASSO能夠在訓練模型時同時選擇因子,并通過L1算法避免過擬合問題[1],最后分析線性回歸和貝葉斯網絡模型在軟件缺陷預測方面的優(yōu)缺點。
在軟件缺陷數(shù)的預測中,因子是離散數(shù)據(jù)和連續(xù)數(shù)據(jù)結合的,采用線性回歸和貝葉斯網絡進行預測,能夠得出這種混合數(shù)據(jù)更適應于用哪種模型來預測。
線性回歸是采用歷史數(shù)據(jù)來預測未來結果的最直接和有效的方法,一般的現(xiàn)實問題都可以用線性回歸模型來描述和預測,且結果往往很準確。
貝葉斯網絡是一個有向無環(huán)圖,圖中的節(jié)點代表一個獨立屬性,節(jié)點的邊代表兩個獨立屬性之間的依賴關系[2]。如圖1所示是一個貝葉斯網絡示例。
圖1 貝葉斯網絡
貝葉斯網絡結構確定時,其中所有屬性的聯(lián)合概率分布作為監(jiān)督值,可以處理有監(jiān)督學習的數(shù)據(jù)仿真及預測;當網絡結構不能事先明確,需要評價函數(shù)來學習得到網絡模型,評價函數(shù)常見的有AIC、BIC以及極大似然估計,都能夠基于假設近似的得到一個較為準確的模型[3]。當網絡結構確立后,就可以通過已知屬性預測未知屬性,例如通過吉布斯采樣采集已知屬性樣本,在網絡中計算未知屬性的后驗概率。
LASSO是一種常用的嵌入式特征選擇方法,用來解決特征值較多,LASSO方法有嵌入式特征選擇方法的優(yōu)點,即能夠同時訓練模型和因子選擇,LASSO相較其它嵌入式特征因子,選擇的另一個顯著優(yōu)勢是通過L1正則化形成了稀疏矩陣[4],能夠最小化的對相關因子進行降維,且訓練結果較為準確,常用于復雜因子模型構建前的數(shù)據(jù)預處理。
缺陷影響因子是在軟件生命周期過程中作為輸入能夠影響輸出結果并引起缺陷的因素。本文所述的因子是從預測對象形成之前的階段中得到的,例如預測需求規(guī)格說明的缺陷數(shù),其缺陷影響因子可能是需求分析人員的水平、用戶的表達程度、需求規(guī)格說明文檔的規(guī)模等,由此可以得知缺陷影響因子的值既有離散數(shù)也有連續(xù)數(shù)。
由于一個軟件缺陷在未確定原因之前,可能有多個影響因子[10],因此借鑒故障樹分析的思路,采用構建無環(huán)有向圖的方式來分析軟件的缺陷影響因子。
因在軟件生命周期過程活動中引入的各種因素具有直接、間接兩種關系,沒有循環(huán)關系,因此可以構建有向無環(huán)圖。
具體采用兩種構建方式,一種是正向構建,即根據(jù)經驗通過原因推導缺陷,一種是逆向構建,通過常用的故障模式來逆向推導原因。構建有向無環(huán)圖默認兩個前提,一是缺陷和故障模式是一個抽象的概念,即在圖中只能作為一個節(jié)點;二是圖中的原因能夠量化,一般對定量的因子直接采集數(shù)據(jù),對定性的因子則采用專家法進行量化。
正向構建時,是由原因推導出可能的結果,構建步驟如下:
從軟件生命周期活動中選取一個活動或一個產品作為起點v0;
分析起點與缺陷相關的能夠量化的因素,每個因素作為有向無環(huán)圖的起點w0i(i=1,2,3…n,n為因素個數(shù));
以v0為輸入,分析其輸出活動或產品,記為v1i(i=1,2,3…m,m為輸出個數(shù));
分析v1i,得到所有可能與缺陷相關且能夠量化的因素,記為w1j(j=1,2,3…l,l為因素個數(shù));
分析w0i與w1j,如果存在相關關系則繪出由w0i到w0j的路徑;
循環(huán)上述步驟構建有向無環(huán)圖V和W,直到目標節(jié)點,即缺陷節(jié)點。
上述步驟中,缺陷節(jié)點的輸入可以是直接影響因子,也可以是間接影響因子,隨應用經驗和分析粒度決定。上述步驟構建的W圖中所有的節(jié)點wij即缺陷影響因子,具體W圖如圖2所示。
圖2 缺陷影響因子關系有向無環(huán)W圖
逆向構建網絡步驟與正向構建思想相同,思路相反,這里不再贅述。為后續(xù)計算方便,圖2用表格的形式表示見表1。
表1 W圖的表格化形式
注:節(jié)點表示因子,前驅表示輸入,后繼表示輸出,邊表示與后繼的相關度。
將上節(jié)所構建的有向無環(huán)圖進行模型化處理即形成貝葉斯網絡模型,進而求解最優(yōu)貝葉斯網絡作為預測模型。
求最優(yōu)貝葉斯網絡是一個NP難題,常用的算法有貪心算法和約束消減算法,因缺陷預測前不確定影響因子的重要性,采用貪心算法結合LASSO特征選擇可以避免重要因子被剔除,貪心算法需要一個初始的網絡模型[11]。
假設有N個缺陷影響因子,構建網絡為B,因子參數(shù)為θ,將網絡B的結構記為
B=
(1)
其中,G為一個有向無環(huán)圖,如圖2所示,每一個節(jié)點對應一個影響因子,參數(shù)θ用來描述因子之間的關系,假設屬性xi在圖G中的父節(jié)點為πi,則將每個屬性的條件概率結構記為
θxi/yi=PB(xi/π)
(2)
假設G有n個節(jié)點,則網絡G的屬性x1,x2,x3…xn的聯(lián)合概率分布定義為
(3)
在式(3)中,πi為第i個類別,訓練集為D,則有
(4)
式中:Dπi,xi為第i個屬性在第i個分類里的個數(shù),Dπi為第i個分類的個數(shù)。依據(jù)經驗僅僅比例估計條件概率結果較為脆弱,尤其當屬性較多訓練數(shù)據(jù)較少的情況下,因此將式(4)改為如下形式
(5)
式(5)采用m估計方法,其中m是等價樣本大小的參數(shù),p為用戶指定的參數(shù),一般依據(jù)m*p=1對p進行指定。
屬性量化值為連續(xù)型數(shù)值,則
(6)
(7)
為避免屬性和類別混淆,式(7)中πj為訓練集D中第j個類別,m為類別個數(shù),n為屬性個數(shù),當預測目標是一個值時,m=1。
貝葉斯網絡一般的學習算法是針對初始網絡模型定義一個評分函數(shù),求評分函數(shù)的最大或最小值進而得到訓練模型[12]。本文確定訓練模型的思路是將初始貝葉斯網絡模型與LASSO特征選擇方法結合后進行學習,在得到一個系數(shù)稀疏解的同時將訓練模型確定下來的過程。
(1)確定模型
按照LASSO特征選擇方法,線性回歸以及概率函數(shù)以平方損失函數(shù)為評估模型,優(yōu)化目標隨屬性的離散和連續(xù)特征不同而不同[11]。依據(jù)經驗,結合上節(jié)的初始網絡模型構建優(yōu)化目標評估函數(shù)如下
(8)
式(8)中,θT為矩陣θ的轉置矩陣,λ為正則化參數(shù),在LASSO中規(guī)定λ>0。
(2)求解模型參數(shù)
求解目標函數(shù)是求解使式(8)的值最小的θ的過程,在計算過程中最小二乘法是數(shù)學求解過程,但需要x滿秩,而在缺陷預測中并不能保證訓練集的每個因子都有量化值,所以采用梯度下降算法進行求解。
(9)
然后帶入梯度下降算法,梯度下降算法如下:
輸入:數(shù)據(jù)集D;
特征集A;
評估函數(shù)f(θ);//式(8)和式(9),y和x作為常量
步長step;//根據(jù)經驗取值,會影響迭代收斂速度
迭代閾值ξ
正則化參數(shù)λ;
過程:
初始化變量
f_change=f(θ);//兩次殘差的差值
f_current=f(θ);//當前殘差
循環(huán)求解
While D 非空 do
//梯度下降θ的值
θ=θ-step*f(θ′);//依據(jù)經驗f(θ′) 為f(θ) 的導數(shù)
//計算殘差差值
f_change=f_current-f(θ);
//計算殘差
f_current=f(θ);
end while
將θ加入到中;
將該步的數(shù)據(jù)對y,x從D中刪除;
End while
(3)缺陷預測
本文采用某單位近五年積累的軟件研制過程的數(shù)據(jù)作為訓練和測試樣本,采集項目類型為辦公自動化項目需求規(guī)格說明的缺陷數(shù),以及其它因子。生產這些數(shù)據(jù)的研發(fā)過程是嚴格按照GJB5000A-2008《軍用軟件研制能力成熟度模型》中的規(guī)定而制定和執(zhí)行的,因此這些數(shù)據(jù)具有普遍性和代表性。具體如圖3所示。
圖3 需求缺陷影響因子集合
由于圖3中的缺陷因子在某單位的數(shù)據(jù)積累中都是定量的數(shù)據(jù),可以直接采集到。
對上述各屬性進行特征值選擇,采用R語言中的glmnet包進行,主要代碼如下:
#lasso特征選擇代碼
mydata<-read.csv(file="C:\Users\user17\Downloads\樣本數(shù)據(jù).csv",header=TRUE)
x=as.matrix(mydata[, 1:8])
y=as.matrix(mydata[, 9])
la<-lars(x,y,type="lar")
plot(la)
summary(la)
運行結果如下:
LAR sequence
Computing X′X…
LARS Step 1: Variable 2 added
LARS Step 2: Variable 3 added
LARS Step 3: Variable 4 added
LARS Step 4: Variable 5 added
LARS Step 5: Variable 6 added
LARS Step 6: Variable 1 added
LARS Step 7: Variable 7 added
LARS Step 8: Variable 8 added
LARS Step 9: Variable 9 added
Computing residuals, RSS etc…
LARS/LAR
Call: lars(x=x,y=y, type="lar", trace=10)
Df Rss Cp
0 1 178.15 30.2213
1 2 165.24 16.2195
2 3 160.34 15.7000
3 4 160.11 13.6305
4 5 158.34 9.0796
5 6 156.87 6.6245
6 7 154.23 4.8069
7 8 148.40 8.0000
7 9 132.40 10.0000
依據(jù)運行結果,選擇CP值最小的步驟,即選擇屬性x2,x3,x4,x5,x6,x1,x7, 去掉x8,x9。
依據(jù)上一節(jié)結果嗎,去掉了x中與y相關性弱的節(jié)點,采用R語言進行線性回歸,代碼如下:
predata<-read.table("E:\mayou\test.csv",header=FALSE,sep=",")
y<-mydata[,8]
x1<-mydata[,1]
x2<-mydata[,2]
x3<-mydata[,3]
x4<-mydata[,4]
x5<-mydata[,5]
x6<-mydata[,6]
x7<-mydata[,7]
plot<-lm(y~x1+x2+x3+x4+x5+x6+x7,data=predata)
plot
#結果
Call:
lm(formula=y~x1+x2+x3+x4+x5+x6+x7,data=predata)
Coefficients:
(Intercept) x1 x2 x3 x4
42.437 42.4607 0.9408 -4.4549 -2.9723
x5 x6 x7
4.1206 -0.7707 -8.2829
#回歸診斷
plot(plot)
從運行結果可知,最后的模型為:
Y=2.4607x1+0.9408x2-4.4549x3-2.9723x4+4.1206x5-0.7707x6-8.2829x7, 依據(jù)上述預測模型,假設置信區(qū)間為0.95,將測試集帶入,代碼和結果如下:
predata<-read.table("E:\mayou\test.csv",header=FALSE,sep=",")
plot<-lm(v8~V1+V2+V3+V4+V5+V6+V7,data=predata)
newdata<-read.table("E:\mayou\test1.csv",header=FALSE,sep=",")
predict(plot,data.frame(newdata[,1:7]))
#結果
1 2 3 4
13.6510211 -0.6531836 13.6781614 26.7401970
5 6 7 8
7.7557199 23.2157111 31.3797127 22.2473188
9 10 11
5.0710729 24.7817223 2.4761034
與真實值的對比結果如下:
323027257131092243793326323518202621
誤差均值為2.8。
根據(jù)各因子的潛在關系,構建貝葉斯網絡模型,如圖4所示。
圖4 貝葉斯網絡模型
采集圖4中所有節(jié)點的數(shù)據(jù),按照3∶1的比例分為訓練數(shù)據(jù)和測試數(shù)據(jù)集合。基于訓練數(shù)據(jù),用R語言繪制貝葉斯網絡模型[12],計算得到因子系數(shù)θ的值,主要代碼如下:
#安裝R貝葉斯庫
install.packages("D:/bnlearn_4.1.zip", repos = NULL, type = "source")
#引入庫
library(bnlearn)
#導入數(shù)據(jù)
mydata1<-read.table("D:/Rworkspace/sample.csv",header=FALSE,sep=",")
#離散化缺陷數(shù)
c8<-t(v8)
cc8<-sort(c8)
group<-c(seq(1,max(cc8),10),max(cc8))
t8<-cut(c8,breaks=group,labels=1:(length(group)-1))
t8
[1] 3 2 1 1 4 1 2 2 2 4 4 2 1 1 2 6 5 1 3 2 1 1 1 1 3 3 4 6 5 4 7 11 1 2 3
[36] 2 5 5 2 1 4 1 1
netdata<-cbind(predata[,1:7],t8)
#繪制貝葉斯節(jié)點
bayesnet<-empty.graph(names(mydata))
bayesnet<-set.arc(bayesnet,"V1","V2")
bayesnet<-set.arc(bayesnet,"V3","V2")
bayesnet<-set.arc(bayesnet,"V1","V3")
bayesnet<-set.arc(bayesnet,"V1","V4")
bayesnet<-set.arc(bayesnet,"V3","V4")
bayesnet<-set.arc(bayesnet,"V6","V4")
bayesnet<-set.arc(bayesnet,"V1","V8")
bayesnet<-set.arc(bayesnet,"V2","V8")
bayesnet<-set.arc(bayesnet,"V3","V8")
bayesnet<-set.arc(bayesnet,"V4","V8")
bayesnet<-set.arc(bayesnet,"V5","V8")
bayesnet<-set.arc(bayesnet,"V6","V8")
bayesnet<-set.arc(bayesnet,"V7","V8")
plot(bayesnet)
fitted<-bn.fit(bayesnet,mydata,method="mle")
fitted
pre<-predict(fitted,data=mydata,node='V8')
#顯示網絡
plot(bayesnet)
#用導入的數(shù)據(jù)訓練模型
fitted<-bn.fit(bayesnet,mydata1,method="mle")
#預測
newdata<-read.table("E:\mayou\test1.csv",header=FALSE,sep=",")
pre<-predict(fitted,data=newdata,node=′V8′)
運行結果中,需求規(guī)格缺陷節(jié)點的概率結果為:
Parameters of node V8 (Gaussian distribution)
Conditional density: V8|V1+V2+V3+V4+V5+V6+V7
Coefficients:
(Intercept) V1 V2
4.47062148 0.25553858 0.09860272
V3 V4 V5
-0.44597170 -0.33725753 0.40965399
V6 V7
-0.02004204 -0.81317686
Standard deviation of the residuals: 1.466314
由于缺陷數(shù)據(jù)是連續(xù)值,進行了離散化,因此需要將預測結果的離散值轉換為缺陷數(shù),本文以最大閾值為預測值,與真實結果對比見表2。
表2 貝葉斯網絡預測結果與真實結果對比
貝葉斯網絡預測結果誤差均值為:3.5。
由上述兩個實例結果可知,兩個模型都比較適合預測缺陷數(shù)量,但由于缺陷數(shù)是連續(xù)值,貝葉斯網絡預測模型的準確度要低一些,且計算較復雜。
由于在構建貝葉斯網絡時具有極大的主觀性,因此當數(shù)據(jù)本身已經是二次加工后的,例如定性數(shù)據(jù)轉換為定量數(shù)據(jù)這種情況,則貝葉斯網絡的預測準確度就會明顯降低。因此,當預測樣本數(shù)據(jù)經過二次加工后形成,且具有主觀性的情況下,采用線性模型預測結果可能更準確一些。
軟件缺陷由于多種原因導致,在實踐中由于各種原因量化的不準確性等問題,理論預測的結果往往不盡如人意,要解決類似的問題,應該從數(shù)據(jù)積累入手,逐步尋找規(guī)律提升預測的準確性。