In this notebook we will compare bnlearn's hillclimbing algorithm and tabu algorithm against Gobnilp. We'll do this with default settings and BIC scoring. First let's use bnlearn's hill-climbing algorithm to learn from its built-in alarm dataset.
library(bnlearn)
data(alarm)
alarmhc.bn <- hc(alarm)
alarmhc.bn
score(alarmhc.bn,alarm,type="bic")
Attaching package: ‘bnlearn’ The following object is masked from ‘package:stats’: sigma
Bayesian network learned via Score-based methods model: [HIST][HRBP][PAP][FIO2][ANES][ERCA][LVF|HIST][PMB|PAP][ERLO|HRBP][PCWP|LVF] [HR|HRBP:ERLO][HREK|HR:ERCA][HRSA|HR:ERCA][LVV|PCWP:LVF][CCHL|HR][CVP|LVV] [MINV|CCHL][STKV|LVF:LVV][CO|STKV:HR][HYP|LVV:STKV][VALV|MINV][INT|MINV:VALV] [PVS|FIO2:VALV][ACO2|CCHL:VALV][PRSS|INT:VALV][SHNT|PMB:INT] [VLNG|MINV:INT:VALV][SAO2|SHNT:PVS][ECO2|ACO2:VLNG][KINK|PRSS:VLNG] [VTUB|PRSS:MINV:INT][TPR|SAO2:CCHL][DISC|VTUB][BP|TPR:CO][APL|TPR] [VMCH|DISC:VTUB][MVS|VMCH] nodes: 37 arcs: 53 undirected arcs: 0 directed arcs: 53 average markov blanket size: 3.46 average neighbourhood size: 2.86 average branching factor: 1.43 learning algorithm: Hill-Climbing score: BIC (disc.) penalization coefficient: 4.951744 tests used in the learning procedure: 2718 optimized: TRUE
OK, so hill-climbing has quickly found a BN with BIC score of -220761.687712649. Let's see how tabu search does.
alarmtabu.bn <- tabu(alarm)
alarmtabu.bn
score(alarmtabu.bn,alarm,type="bic")
Bayesian network learned via Score-based methods model: [PAP][FIO2][LVF][ANES][ERLO][HR][ERCA][HIST|LVF][HRBP|ERLO:HR][HREK|HR:ERCA] [HRSA|HR:ERCA][PMB|PAP][LVV|LVF][CCHL|HR][CVP|LVV][PCWP|LVV][MINV|CCHL] [STKV|LVF:LVV][CO|STKV:HR][HYP|LVV:STKV][VALV|MINV][INT|MINV:VALV] [PVS|FIO2:VALV][ACO2|CCHL:VALV][PRSS|INT:VALV][SHNT|PMB:INT] [VLNG|MINV:INT:VALV][SAO2|SHNT:PVS][ECO2|ACO2:VLNG][KINK|PRSS:VLNG] [VTUB|PRSS:MINV:INT][TPR|SAO2:CCHL][DISC|VTUB][BP|TPR:CO][APL|TPR] [VMCH|DISC:VTUB][MVS|VMCH] nodes: 37 arcs: 51 undirected arcs: 0 directed arcs: 51 average markov blanket size: 3.41 average neighbourhood size: 2.76 average branching factor: 1.38 learning algorithm: Tabu Search score: BIC (disc.) penalization coefficient: 4.951744 tests used in the learning procedure: 3759 optimized: TRUE
Tabu search is also very quick and finds a slightly better network with a score of -220727.331509064.
Now let's see how Gobnilp does on the same data.
library(reticulate)
m <- import("pygobnilp.gobnilp")$Gobnilp()
m$learn(alarm,plot=FALSE,score='DiscreteBIC')
m$learned_bn
********** BN has score -218632.31567964845 ********** ACO2<-CCHL,VALV -7365.911141511497 VALV<-INT,VLNG -3656.8027442121847 CCHL<-SAO2,TPR -5079.249250888837 ANES<- -6708.582939203498 APL<- -1134.1626027276418 BP<-CO,TPR -10757.751590048356 TPR<-APL -21601.664138686014 CO<-HR,STKV -5340.866593873159 SAO2<-PVS,SHNT -2469.044793071091 HR<-CCHL -7361.101683937812 STKV<-HYP,LVF -9121.2931288759 CVP<-LVV -6151.495470079891 LVV<-HYP,LVF -7523.033824999751 DISC<- -6447.082832428477 ECO2<-ACO2,VLNG -3540.04821713251 VLNG<-INT,VTUB -7559.536030079053 ERCA<- -6480.204445349012 ERLO<- -3981.1433847172025 FIO2<- -4007.581975260592 HIST<-LVF -1488.3560086436557 LVF<- -3945.7594222796324 HRBP<-ERLO,HR -2822.433256557475 HREK<-ERCA,HR -3136.739214776165 HRSA<-ERCA,HR -3236.139444268571 HYP<-APL -10075.035674773593 INT<- -6658.474123732305 KINK<-VLNG,VTUB -2857.4047588872545 VTUB<-DISC,VMCH -3418.0624644648706 MINV<-INT,VLNG -3909.141376477159 MVS<- -8115.930179708217 PAP<- -8357.309503322173 PCWP<-LVV -4417.152413059788 PMB<-PAP -782.1722102996109 PRSS<-INT,KINK,VTUB -17104.845843232193 PVS<-FIO2,VALV -1881.4613201566744 SHNT<-INT,PMB -3967.786877754141 VMCH<-MVS -6171.55480017248 ********** bnlearn modelstring = [ACO2|VALV:CCHL][VALV|INT:VLNG][CCHL|TPR:SAO2][ANES][APL][BP|TPR:CO][TPR|APL][CO|HR:STKV][SAO2|SHNT:PVS][HR|CCHL][STKV|LVF:HYP][CVP|LVV][LVV|LVF:HYP][DISC][ECO2|VLNG:ACO2][VLNG|INT:VTUB][ERCA][ERLO][FIO2][HIST|LVF][LVF][HRBP|ERLO:HR][HREK|ERCA:HR][HRSA|ERCA:HR][HYP|APL][INT][KINK|VLNG:VTUB][VTUB|DISC:VMCH][MINV|INT:VLNG][MVS][PAP][PCWP|LVV][PMB|PAP][PRSS|INT:VTUB:KINK][PVS|VALV:FIO2][SHNT|INT:PMB][VMCH|MVS]
Gobnilp is a lot slower, but on the other hand, the learned BN has a much better BIC score of -218632.30999687716. Note that Gobnilp is being run here with its default limit of at most 3 parents for any node. It could be that an even higher-scoring network exists where some nodes have more than 3 parents.
We can get bnlearn to score the network learned by Gobnilp to check that we get the same BIC score.
score(model2network(m$learned_bn$bnlearn_modelstring()),alarm,type="bic")
Sure enough we get the same score (modulo some numerical imprecision). Now let's run the same experiment with bnlearn's built-in hailfinder dataset
data(hailfinder)
hailfinderhc.bn <- hc(hailfinder)
hailfinderhc.bn
score(hailfinderhc.bn,hailfinder,type="bic")
hailfindertabu.bn <- tabu(hailfinder)
hailfindertabu.bn
score(hailfindertabu.bn,hailfinder,type="bic")
m$learn(hailfinder,plot=FALSE,score='DiscreteBIC')
m$learned_bn
score(model2network(m$learned_bn$bnlearn_modelstring()),hailfinder,type="bic")
Bayesian network learned via Score-based methods model: [N07muVerMo][SubjVertMo][QGVertMotion][SatContMoist][RaoContMoist] [VISCloudCov][IRCloudCover][AMInstabMt][WndHodograph][MorningBound] [LoLevMoistAd][Date][MorningCIN][LIfr12ZDENSd][AMDewptCalPl][LatestCIN][LLIW] [CombVerMo|N07muVerMo:SubjVertMo:QGVertMotion] [CombMoisture|SatContMoist:RaoContMoist][CombClouds|VISCloudCov:IRCloudCover] [Scenario|Date][CurPropConv|LatestCIN:LLIW][AreaMesoALS|CombVerMo] [AreaMoDryAir|CombVerMo:CombMoisture][ScenRelAMCIN|Scenario] [ScenRelAMIns|Scenario][ScenRel34|Scenario][ScnRelPlFcst|Scenario] [Dewpoints|Scenario][LowLLapse|Scenario][MeanRH|Scenario][MidLLapse|Scenario] [MvmtFeatures|Scenario][RHRatio|Scenario][SfcWndShfDis|Scenario] [SynForcng|Scenario][TempDis|Scenario][WindAloft|Scenario] [WindFieldMt|Scenario][WindFieldPln|Scenario] [CldShadeOth|CombVerMo:AreaMoDryAir:CombClouds] [AMCINInScen|ScenRelAMCIN:MorningCIN] [AMInsWliScen|ScenRelAMIns:LIfr12ZDENSd:AMDewptCalPl] [InsInMt|CldShadeOth:AMInstabMt][OutflowFrMt|InsInMt:WndHodograph] [CldShadeConv|InsInMt:WndHodograph][MountainFcst|InsInMt] [Boundaries|WndHodograph:OutflowFrMt:MorningBound] [CompPlFcst|CombVerMo:CldShadeOth:CldShadeConv][CapChange|CompPlFcst] [InsChange|CompPlFcst:LoLevMoistAd][CapInScen|CompPlFcst:AMCINInScen] [InsSclInScen|InsChange:AMInsWliScen] [PlainsFcst|ScenRelAMCIN:InsSclInScen:CurPropConv] [N34StarFcst|ScenRel34:PlainsFcst][R5Fcst|MountainFcst:N34StarFcst] nodes: 56 arcs: 64 undirected arcs: 0 directed arcs: 64 average markov blanket size: 3.25 average neighbourhood size: 2.29 average branching factor: 1.14 learning algorithm: Hill-Climbing score: BIC (disc.) penalization coefficient: 4.951744 tests used in the learning procedure: 5060 optimized: TRUE
Bayesian network learned via Score-based methods model: [N07muVerMo][SubjVertMo][QGVertMotion][SatContMoist][RaoContMoist] [VISCloudCov][IRCloudCover][AMInstabMt][WndHodograph][MorningBound] [LoLevMoistAd][Date][MorningCIN][LIfr12ZDENSd][AMDewptCalPl][LatestCIN][LLIW] [CombVerMo|N07muVerMo:SubjVertMo:QGVertMotion] [CombMoisture|SatContMoist:RaoContMoist][CombClouds|VISCloudCov:IRCloudCover] [Scenario|Date][CurPropConv|LatestCIN:LLIW][AreaMesoALS|CombVerMo] [AreaMoDryAir|CombVerMo:CombMoisture][ScenRelAMCIN|Scenario] [ScenRelAMIns|Scenario][ScenRel34|Scenario][ScnRelPlFcst|Scenario] [Dewpoints|Scenario][LowLLapse|Scenario][MeanRH|Scenario][MidLLapse|Scenario] [MvmtFeatures|Scenario][RHRatio|Scenario][SfcWndShfDis|Scenario] [SynForcng|Scenario][TempDis|Scenario][WindAloft|Scenario] [WindFieldMt|Scenario][WindFieldPln|Scenario] [CldShadeOth|CombVerMo:AreaMoDryAir:CombClouds] [AMCINInScen|ScenRelAMCIN:MorningCIN] [AMInsWliScen|ScenRelAMIns:LIfr12ZDENSd:AMDewptCalPl] [InsInMt|CldShadeOth:AMInstabMt][OutflowFrMt|InsInMt:WndHodograph] [CldShadeConv|InsInMt:WndHodograph][MountainFcst|InsInMt] [Boundaries|WndHodograph:OutflowFrMt:MorningBound] [CompPlFcst|CombVerMo:CldShadeOth:CldShadeConv][CapChange|CompPlFcst] [InsChange|CompPlFcst:LoLevMoistAd][CapInScen|CompPlFcst:AMCINInScen] [InsSclInScen|InsChange:AMInsWliScen] [PlainsFcst|ScenRelAMCIN:InsSclInScen:CurPropConv] [N34StarFcst|ScenRel34:PlainsFcst][R5Fcst|MountainFcst:N34StarFcst] nodes: 56 arcs: 64 undirected arcs: 0 directed arcs: 64 average markov blanket size: 3.25 average neighbourhood size: 2.29 average branching factor: 1.14 learning algorithm: Tabu Search score: BIC (disc.) penalization coefficient: 4.951744 tests used in the learning procedure: 6014 optimized: TRUE
********** BN has score -990162.8050400699 ********** AMCINInScen<-MorningCIN,ScenRelAMCIN -14827.804016810527 ScenRelAMCIN<-Scenario -54.46918153893982 MorningCIN<- -22451.66219109244 AMDewptCalPl<- -21401.356146585407 AMInsWliScen<-AMDewptCalPl,ScenRelAMIns -19662.46130603437 ScenRelAMIns<-ScnRelPlFcst -272.34590769473465 AMInstabMt<- -21979.829675807607 AreaMesoALS<-CombVerMo -59.42092531521676 CombVerMo<-N07muVerMo,QGVertMotion,SubjVertMo -16140.95428595368 AreaMoDryAir<-CombMoisture,CombVerMo -15243.049284360612 CombMoisture<-RaoContMoist,SatContMoist -22201.123014157503 Boundaries<-MorningBound,OutflowFrMt,WndHodograph -15464.935101862531 OutflowFrMt<-InsInMt,WndHodograph -12684.806287718633 WndHodograph<- -27534.501610060404 MorningBound<- -19880.990986381013 CapChange<-CompPlFcst -29.71046265760838 CompPlFcst<-CldShadeConv,CldShadeOth,CombVerMo -20195.300816671996 CapInScen<-AMCINInScen,CapChange -8769.900372311302 CldShadeConv<-InsInMt,WndHodograph -11930.685819990145 InsInMt<-AMInstabMt,CldShadeOth -9186.457382838727 CldShadeOth<-AreaMoDryAir,CombClouds,CombVerMo -10808.429735453852 CombClouds<-IRCloudCover,VISCloudCov -12416.50118390726 VISCloudCov<- -18807.3986119801 IRCloudCover<- -20081.087288850224 RaoContMoist<- -26451.36187705046 SatContMoist<- -26381.105022929798 QGVertMotion<- -24641.079933873705 SubjVertMo<- -24741.633744298222 N07muVerMo<- -27740.002999058954 CurPropConv<-PlainsFcst -26471.80989062072 PlainsFcst<-InsSclInScen,Scenario -17466.172736436805 Date<-Scenario -33427.328674309276 Scenario<-ScnRelPlFcst -544.6918153894782 Dewpoints<-Scenario -30495.957816711856 InsChange<-CompPlFcst,LoLevMoistAd -15438.63271567674 LoLevMoistAd<- -26697.017932815765 InsSclInScen<-AMInsWliScen,InsChange -11294.899810438172 LIfr12ZDENSd<-AMDewptCalPl,AMInsWliScen -21266.227971706914 LLIW<-CurPropConv -21747.034006031623 LatestCIN<-CurPropConv,LLIW -19532.950433050428 LowLLapse<-ScnRelPlFcst -22799.692963671845 ScnRelPlFcst<-WindFieldMt -43976.15779083476 MeanRH<-Scenario -16858.5681959753 MidLLapse<-Scenario -18464.904348591688 MountainFcst<-InsInMt -15270.69659773143 MvmtFeatures<-Scenario -19160.848763138114 N34StarFcst<-PlainsFcst,ScenRel34 -5277.188589205412 ScenRel34<-ScnRelPlFcst -217.8767261557948 R5Fcst<-MountainFcst,N34StarFcst -89.1313879728207 RHRatio<-ScnRelPlFcst -16934.816498959048 WindFieldMt<- -13845.133818871916 SfcWndShfDis<-ScnRelPlFcst -28237.61080138319 SynForcng<-ScnRelPlFcst -27008.48382226494 TempDis<-Scenario -20415.645867981773 WindAloft<-ScnRelPlFcst -19892.419524674344 WindFieldPln<-Scenario -25290.54036622376 ********** bnlearn modelstring = [AMCINInScen|ScenRelAMCIN:MorningCIN][ScenRelAMCIN|Scenario][MorningCIN][AMDewptCalPl][AMInsWliScen|ScenRelAMIns:AMDewptCalPl][ScenRelAMIns|ScnRelPlFcst][AMInstabMt][AreaMesoALS|CombVerMo][CombVerMo|QGVertMotion:SubjVertMo:N07muVerMo][AreaMoDryAir|CombVerMo:CombMoisture][CombMoisture|RaoContMoist:SatContMoist][Boundaries|OutflowFrMt:WndHodograph:MorningBound][OutflowFrMt|WndHodograph:InsInMt][WndHodograph][MorningBound][CapChange|CompPlFcst][CompPlFcst|CldShadeConv:CombVerMo:CldShadeOth][CapInScen|CapChange:AMCINInScen][CldShadeConv|WndHodograph:InsInMt][InsInMt|AMInstabMt:CldShadeOth][CldShadeOth|AreaMoDryAir:CombClouds:CombVerMo][CombClouds|VISCloudCov:IRCloudCover][VISCloudCov][IRCloudCover][RaoContMoist][SatContMoist][QGVertMotion][SubjVertMo][N07muVerMo][CurPropConv|PlainsFcst][PlainsFcst|Scenario:InsSclInScen][Date|Scenario][Scenario|ScnRelPlFcst][Dewpoints|Scenario][InsChange|LoLevMoistAd:CompPlFcst][LoLevMoistAd][InsSclInScen|InsChange:AMInsWliScen][LIfr12ZDENSd|AMDewptCalPl:AMInsWliScen][LLIW|CurPropConv][LatestCIN|CurPropConv:LLIW][LowLLapse|ScnRelPlFcst][ScnRelPlFcst|WindFieldMt][MeanRH|Scenario][MidLLapse|Scenario][MountainFcst|InsInMt][MvmtFeatures|Scenario][N34StarFcst|ScenRel34:PlainsFcst][ScenRel34|ScnRelPlFcst][R5Fcst|N34StarFcst:MountainFcst][RHRatio|ScnRelPlFcst][WindFieldMt][SfcWndShfDis|ScnRelPlFcst][SynForcng|ScnRelPlFcst][TempDis|Scenario][WindAloft|ScnRelPlFcst][WindFieldPln|Scenario]
Gobnilp takes a very long time to find a BN mainly due to this Python version of Gobnilp computing the necessary local scores slowly. The found BN has a BIC score of -990162.80504007 in contrast to -990474.753096428 BIC score achieved by hill-climbing and tabu.
Finally, let's have a look at the Insurance dataset
data(insurance)
insurancehc.bn <- hc(insurance)
insurancehc.bn
score(insurancehc.bn,insurance,type="bic")
insurancetabu.bn <- tabu(insurance)
insurancetabu.bn
score(insurancetabu.bn,insurance,type="bic")
m$learn(insurance,plot=FALSE,score='DiscreteBIC')
m$learned_bn
score(model2network(m$learned_bn$bnlearn_modelstring()),insurance,type="bic")
Bayesian network learned via Score-based methods model: [RuggedAuto][MakeModel|RuggedAuto][CarValue|RuggedAuto:MakeModel] [Mileage|CarValue][VehicleYear|MakeModel:Mileage:CarValue] [SocioEcon|VehicleYear:MakeModel][Antilock|VehicleYear:MakeModel] [Airbag|VehicleYear:MakeModel][ThisCarDam|Mileage:Antilock] [OtherCar|SocioEcon][Cushioning|RuggedAuto:Airbag] [Accident|ThisCarDam:RuggedAuto][ThisCarCost|ThisCarDam:CarValue] [MedCost|ThisCarDam:Cushioning][DrivQuality|Accident:Mileage] [Theft|ThisCarDam:ThisCarCost][OtherCarCost|RuggedAuto:Accident] [ILiCost|Accident][Age|SocioEcon:DrivQuality] [PropCost|ThisCarCost:OtherCarCost][GoodStudent|Age:SocioEcon] [SeniorTrain|Age:DrivQuality][RiskAversion|Age:DrivQuality:SeniorTrain] [DrivingSkill|RiskAversion:DrivQuality][HomeBase|SocioEcon:RiskAversion] [AntiTheft|SocioEcon:RiskAversion][DrivHist|RiskAversion:DrivingSkill] nodes: 27 arcs: 50 undirected arcs: 0 directed arcs: 50 average markov blanket size: 4.44 average neighbourhood size: 3.70 average branching factor: 1.85 learning algorithm: Hill-Climbing score: BIC (disc.) penalization coefficient: 4.951744 tests used in the learning procedure: 1755 optimized: TRUE
Bayesian network learned via Score-based methods model: [MakeModel][CarValue|MakeModel][Mileage|CarValue] [VehicleYear|MakeModel:Mileage:CarValue][SocioEcon|VehicleYear:MakeModel] [RuggedAuto|VehicleYear:MakeModel][Antilock|VehicleYear:MakeModel] [Airbag|VehicleYear:MakeModel][ThisCarDam|Mileage:Antilock] [OtherCar|SocioEcon][Cushioning|RuggedAuto:Airbag] [Accident|ThisCarDam:RuggedAuto][ThisCarCost|ThisCarDam:CarValue] [MedCost|ThisCarDam:Cushioning][DrivQuality|Accident:Mileage] [Theft|ThisCarDam:ThisCarCost][OtherCarCost|RuggedAuto:Accident] [ILiCost|Accident][Age|SocioEcon:DrivQuality] [PropCost|ThisCarCost:OtherCarCost][GoodStudent|Age:SocioEcon] [SeniorTrain|Age:DrivQuality][RiskAversion|Age:DrivQuality:SeniorTrain] [DrivingSkill|RiskAversion:DrivQuality][HomeBase|SocioEcon:RiskAversion] [AntiTheft|SocioEcon:RiskAversion][DrivHist|RiskAversion:DrivingSkill] nodes: 27 arcs: 50 undirected arcs: 0 directed arcs: 50 average markov blanket size: 4.44 average neighbourhood size: 3.70 average branching factor: 1.85 learning algorithm: Tabu Search score: BIC (disc.) penalization coefficient: 4.951744 tests used in the learning procedure: 2405 optimized: TRUE
********** BN has score -264708.85975077713 ********** Accident<-Antilock,DrivQuality -12310.195471829667 DrivQuality<-DrivingSkill,RiskAversion -4890.4371039384205 Antilock<-MakeModel,VehicleYear -3057.5537899833384 Age<-RiskAversion,SocioEcon -17309.67690764885 RiskAversion<- -22128.556172909437 SocioEcon<-RiskAversion -21914.967878236253 Airbag<-MakeModel,VehicleYear -4426.614280155772 VehicleYear<-SocioEcon -10055.937028683162 MakeModel<-SocioEcon -18409.07715172664 AntiTheft<-RiskAversion,SocioEcon -5427.903744351533 CarValue<-MakeModel,Mileage,VehicleYear -11811.232170180529 Mileage<-Accident,DrivQuality -23647.681978155444 Cushioning<-Airbag,RuggedAuto -14785.848505297057 RuggedAuto<-MakeModel,VehicleYear -12165.511166709075 DrivHist<-DrivingSkill,RiskAversion -9560.920896777798 DrivingSkill<-Age,SeniorTrain -17818.740523142405 SeniorTrain<-Age,RiskAversion -1712.9733331992413 GoodStudent<-Age,SocioEcon -1937.6149662237774 HomeBase<-RiskAversion,SocioEcon -14317.997554087879 ILiCost<-Accident -2459.8973674972813 MedCost<-Accident,Cushioning -3741.377927855121 OtherCar<-SocioEcon -10846.63435568069 OtherCarCost<-Accident,RuggedAuto -4190.349482680152 PropCost<-OtherCarCost,ThisCarCost -11581.005875670833 ThisCarCost<-CarValue,ThisCarDam -1499.5192944884964 Theft<-ThisCarCost,ThisCarDam -120.29650307184178 ThisCarDam<-Accident,RuggedAuto -2580.3383205964406 ********** bnlearn modelstring = [Accident|DrivQuality:Antilock][DrivQuality|RiskAversion:DrivingSkill][Antilock|VehicleYear:MakeModel][Age|RiskAversion:SocioEcon][RiskAversion][SocioEcon|RiskAversion][Airbag|VehicleYear:MakeModel][VehicleYear|SocioEcon][MakeModel|SocioEcon][AntiTheft|RiskAversion:SocioEcon][CarValue|VehicleYear:MakeModel:Mileage][Mileage|DrivQuality:Accident][Cushioning|RuggedAuto:Airbag][RuggedAuto|VehicleYear:MakeModel][DrivHist|RiskAversion:DrivingSkill][DrivingSkill|SeniorTrain:Age][SeniorTrain|RiskAversion:Age][GoodStudent|Age:SocioEcon][HomeBase|RiskAversion:SocioEcon][ILiCost|Accident][MedCost|Cushioning:Accident][OtherCar|SocioEcon][OtherCarCost|Accident:RuggedAuto][PropCost|OtherCarCost:ThisCarCost][ThisCarCost|ThisCarDam:CarValue][Theft|ThisCarCost:ThisCarDam][ThisCarDam|Accident:RuggedAuto]
So hill-climbing finds a BN with BIC score -266113.028472967, tabu search manages -265500.923940094 and Gobnilp achieves -264708.85975077713. In this case Gobnilp was not too slow either (since Insurance only has 27 nodes and the default limit of at most 3 parents per node is in use.)