ComputersandStructuresxxx(2011)xxx–xxx ContentslistsavailableatScienceDirect Computers and Structures journal homepage: www.elsevier.com/locate/compstruc Bayesian identification of a cracked plate using a population-based Markov Chain Monte Carlo method J.M. Nicholsa,⇑, E.Z. Mooreb, K.D. Murphyb aUSNavalResearchLaboratory,4555OverlookAve.,Washington,DC20375,UnitedStates bUniversityofConnecticut,Dept.ofMechanicalEngineering,Storrs,CT06269-3139,UnitedStates a r t i c l e i n f o a b s t r a c t Articlehistory: Estimatingdamageinstructuralsystemsisachallengingproblemduetothecomplexityofthelikelihood Received21September2010 functiondescribingtheobserveddata.FromaBayesianperspectiveacomplicatedlikelihoodmeanseffi- Accepted25March2011 cientsamplingoftheposteriordistributionisdifficultandstandardMarkovChainMonteCarlosamplers Availableonlinexxxx maynolongerbesufficient.Thisworkdescribesapopulation-basedMarkovChainMonteCarloapproach forefficientsamplingofthedamageparameterposteriordistributions.Theapproachisshowntoaccu- Keywords: ratelyestimatethestateofdamageinacrackedplatestructureusingsimulated,free-decayresponse Bayesianinference data.Theuseofthisapproachinidentifyingstructuraldamagehasnotpreviouslybeenexplored. Population-basedMarkovChainMonte PublishedbyElsevierLtd. Carlo Systemidentification Damageidentification 1.Introduction researchers and solutions proposed. Horibe, for example, used a geneticalgorithm(GA)tofindtheMLEofastructure’sparameters Theproblemofidentifyingdamageinstructuralsystemscanbe inasimplevibrationproblem[1].Likewise,Stulletal.considereda effectivelycastasanestimationproblem.Thatistosay:givenase- GA in generating MLEs of structural parameters in a static shell quenceofobserveddatacollectedfromoneormoresensors,esti- bucklingproblem[2].InthecontextofdamagedetectionPanigrahi matethedamagepresence,location,andextent.Thisinformation et al. used a GA to perform parameter estimation by minimizing isnecessaryforpredictingthefuturestateofthestructureinques- differencesbetweenpredictedandobservedmodalproperties[3] tion,i.e.forprognostics. usingalikelihoodbasedonmodalproperties.Additionalworkby Therearetwomainclassesofestimationapproachesavailable: Hwangetal.alsoemployedacostfunctionbasedonmodalprop- themethodofmaximumlikelihoodandBayesianestimation.The ertiesandusedaGAtoestimatestiffnessinacompositestructure former produces a single ‘‘best’’ estimate, defined as that which [4].Theauthorspointoutthatforcertainproblemsmaximizingthe maximizesthe probability ofhavingobserveda givenset ofdata likelihood does not necessarily require a GA. For example, if en- (e.g.structures’vibrationalresponse).Thelatterviewstheparam- oughisknownabouttheparameterdistributionsaprioriandone etertobeestimatedasarandomvariableandthegoalistocom- has a good initial guess, the perturbation approach of Fonseca bine the likelihood with prior information to produce the etal.[5]orXuetal.[6]canbeusedtoobtainMLEsofastructure’s estimatedprobabilitydensityfunction(PDF)associatedwitheach parameters.Intherelatedfieldofstructuralreliability,researchers parameter. havealsolongrecognizedthechallengesinstructuraloptimization The difficulty with either approach in structural estimation problems. The early work of Kiureghian and Dakessian [7], for problemsisthatthelikelihoodfunction(thecoreofbothestima- example,developedamethodforavoidinglocaloptimainestimat- tion approaches) is often an extremely complicated function of ing failure probability in structures. More recently Guo et al. [8] the structure’s parameters with many maxima. This severely developed an algorithm for efficiently obtaining globally optimal restrictstheuseofstandardoptimizationalgorithms.Forexample, solutionsforstructuraldesignproblemspossessingmultiplelocal asimplegradientascentalgorithmcannotbeusedtofindthemax- optima. imumlikelihoodestimate(MLE)inmanystructuralsystemidenti- Inthispaper,however,theauthorsareinterestedinobtaining fication problems. This has been recognized by a number of the entire parameter probability distribution. This distribution can be used to extract both the parameter estimate and credible ⇑ intervals whichprovide a measure of confidence in the estimate. Correspondingauthor. The authors view well-defined credible intervals as essential in E-mailaddress:[email protected](J.M.Nichols). 0045-7949/$-seefrontmatterPublishedbyElsevierLtd. doi:10.1016/j.compstruc.2011.03.013 Pleasecitethisarticleinpressas:NicholsJMetal.Bayesianidentificationofacrackedplateusingapopulation-basedMarkovChainMonteCarlomethod. ComputStruct(2011),doi:10.1016/j.compstruc.2011.03.013 Report Documentation Page Form Approved OMB No. 0704-0188 Public reporting burden for the collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources, gathering and maintaining the data needed, and completing and reviewing the collection of information. Send comments regarding this burden estimate or any other aspect of this collection of information, including suggestions for reducing this burden, to Washington Headquarters Services, Directorate for Information Operations and Reports, 1215 Jefferson Davis Highway, Suite 1204, Arlington VA 22202-4302. Respondents should be aware that notwithstanding any other provision of law, no person shall be subject to a penalty for failing to comply with a collection of information if it does not display a currently valid OMB control number. 1. REPORT DATE 3. DATES COVERED SEP 2010 2. REPORT TYPE 00-00-2010 to 00-00-2010 4. TITLE AND SUBTITLE 5a. CONTRACT NUMBER Bayesian identification of a cracked plate using a population-based 5b. GRANT NUMBER Markov Chain Monte Carlo method 5c. PROGRAM ELEMENT NUMBER 6. AUTHOR(S) 5d. PROJECT NUMBER 5e. TASK NUMBER 5f. WORK UNIT NUMBER 7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES) 8. PERFORMING ORGANIZATION US Naval Research Laboratory,4555 Overlook REPORT NUMBER Ave,Washington,DC,20375 9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES) 10. SPONSOR/MONITOR’S ACRONYM(S) 11. SPONSOR/MONITOR’S REPORT NUMBER(S) 12. DISTRIBUTION/AVAILABILITY STATEMENT Approved for public release; distribution unlimited 13. SUPPLEMENTARY NOTES 14. ABSTRACT Estimating damage in structural systems is a challenging problem due to the complexity of the likelihood function describing the observed data. From a Bayesian perspective a complicated likelihood means efficient sampling of the posterior distribution is difficult and standard Markov Chain Monte Carlo samplers may no longer be sufficient. This work describes a population-based Markov Chain Monte Carlo approach for efficient sampling of the damage parameter posterior distributions. The approach is shown to accurately estimate the state of damage in a cracked plate structure using simulated, free-decay response data. The use of this approach in identifying structural damage has not previously been explored. 15. SUBJECT TERMS 16. SECURITY CLASSIFICATION OF: 17. LIMITATION OF 18. NUMBER 19a. NAME OF ABSTRACT OF PAGES RESPONSIBLE PERSON a. REPORT b. ABSTRACT c. THIS PAGE Same as 10 unclassified unclassified unclassified Report (SAR) Standard Form 298 (Rev. 8-98) Prescribed by ANSI Std Z39-18 2 J.M.Nicholsetal./ComputersandStructuresxxx(2011)xxx–xxx damage detection application. Optimal decisions regarding the s¼wðhÞþn ð1Þ maintenance and usage of a structure are predicated on where the noise vector n(cid:2)n, i=1,...,T is comprised of T, zero- quantifying the uncertainty in the damage estimate. To this end i mean,iidentrieswithjointprobabilitydistribution the Bayesian estimation framework is conceptually the most appealing.InBayesianestimationonecombinespriorinformation PT 1 (cid:3) ni=2r2 wesittihmathteeAlikneulmihbooerdotforepsreoadruccheertshheadveeshiraeddspuocscteesrsioursipnagraamMeaterr- pHðnÞ¼ð2pr2ÞT=2e i ð2Þ kov-Chain-Monte-Carlo(MCMC)approachtoexploretheposterior Thus,thenoiseisassumedtobestationary,Gaussianwithvariance distributioninstructuralparameterestimationproblems.Numer- r2. The likelihood function is then formed by simply substituting ousworksbyBecketal.(seee.g.[9–11])haveusedMCMCmethods n =s (cid:3)w(h)giving i i i toestimatestructuralparameters.AdditionalworkbyGlaseretal. [12] illustrated the approach in detecting stiffness reduction in 1 (cid:3)PT ðsi(cid:3)wiðhÞÞ2=2r2 beams using static measurements. Zhang and Cho [13] also used pHðnÞ(cid:2)pHðsjhÞ¼ð2pr2ÞT=2e i ð3Þ theMCMCapproachtohelpdesignanevolutionaryalgorithmfor performingsystemidentification. whichistheprobabilityofhavingobservedthesequenceofdatas TheMCMCapproach(tobedescribed)worksbyrepeatedlyper- given the model, defined by parameters h. The aforementioned turbingtheparameter(s)ofinterestandeitheracceptingorreject- MLEs are obtained by maximizing Eq. (16). However, the goal is ing the proposals based on a well-defined criteria. However, the toobtainestimatesofboththeparameterandtheamountofuncer- MCMC algorithm can become ‘‘stuck’’ in a local minima, again taintyintheparameterestimate.ForthisreasontheBayesianesti- due to the complicated likelihood function (see Section 2). This mationphilosophyisadoptedthusprovidingestimatesoftheentire is,ofcourse,whysomanyresearchersuseGAstoexplorecompli- parameter posteriordistributionfromwhichcredibleintervalsfor cated parameter spaces in obtaining MLEs. As a result, recent ourparameterscaneasilybeobtained. workshaveproposedtofuseGAswithMCMCinordertomoreeffi- Bayes’ rule states that the joint posterior distribution of our ciently explore the posterior distribution. Efforts to this end in- modelparametersmaybefoundbytherelationship clude the work of Vrugt et al. [14] who combined differential evolution (a genetic algorithm search technique) with MCMC to p ðhÞ¼pHðsjhÞppðhÞ ð4Þ drawsamplesfromparametersinacomplexsoilmoisturemodel. H pDðsÞ AlsonotedistheworkofZhangandCho[13]whocombinedthe wherepp(h)isthejointpriorparameterdistribution,reflectingany searching capabilities of an evolutionary algorithm with MCMC aprioriinformationonemighthaveaboutourparameter,andp (s) D to produce faster convergence of parameter estimates using data is the joint distribution of our acquired data. Of course the above generatedfromalasersystem.Theuseofmultiplesolutionsthat expression is for the joint parameter distribution whereas one is canexchangeinformationastheyexploretheposteriorhavecome typically interested in the marginal (individual) posteriors, i.e. tobeknownas‘‘population-basedMCMC’’methods.Aniceover- p ðh Þ. Analytically this would require integrating Eq. (4) over viewofpopulation-basedMCMC(Pop-MCMC)methodsisprovided Hp p theotherP(cid:3)1parameters,e.g. byJasraetal.[15]wheretheinitialideawascreditedtoGeyer[16]. Z OtherrecentPop-MCMCworksincludetheimagingapplicationof p ðh Þ¼ p ðhÞdh ð5Þ Kimetal.[17].Althoughthisapproachhasseenrecentuseinthese Hp p RP(cid:3)1 H (cid:3)p and other complicated estimation problems it has not yet been R where the notation dh denotes the multi-dimensional inte- usedinstructuraldynamicsproblemsdespiteobviousadvantages RP(cid:3)1 (cid:3)p graloverallparametersotherthanh .Thiscannottypicallybedone over conventional MCMC in sampling complicated posterior p inclosedformandinsteadtheauthorsresorttonumericalmethods. distributions. The Markov Chain Monte Carlo (MCMC) algorithm was pro- OurgoalinthismanuscriptisthereforetoadoptthePop-MCMC posed by Hastings [18] as a means of drawing samples from Eq. approachtoadifficultparameterestimationprobleminstructural (4)directly(withouthavingtointegrate).TheMCMCalgorithmis dynamics, namely the estimation of parameters governing crack nowfairlystandardandthereaderisreferredto[19]foradescrip- damage in an aluminum plate. The specific implementation of tioninthecontextofstructuraldynamics.Thealgorithmbuildsa thepop-MCMCalgorithmisnewasistheapplicationtostructural Markovchainforeachparametersuchthatsamplesfromthechain vibration problems. The population-based approach is shown to are, in fact, samples from the desired posterior distribution. The providereliableestimatesofthecrackparametersunderrealistic chainisformedbyfirstproposingstatetransitionsforeachparam- levelsofnoise.Theseestimatesarenotalwayspossibleusingstan- eterviaaso-calledproposaldistribution,qðh(cid:4)jh Þ.Thisdistribution dardMCMCaswillbeshown.Additionally,implementationofthe providesaruleforgeneratingacandidateparpampetervalueh(cid:4)given Pop-MCMCalgorithmrequiresthattheforwardmodelbecompu- p thecurrentvalueinthechainh .Theproposedparameter(keeping tationallyefficient.Section4isthereforedevotedtodevelopinga p allotherparametersfixedattheircurrentvalue)isthenaccepted low-dimensionalfiniteelement(FE)crackmodel.Themodeluses orrejectedwithprobability specially tailored elements, resulting in a model that retains the fidelity of the standard FE approach without the computational r ¼min(pHðsjh(cid:4)pÞpppðh(cid:4)pÞ;1): ð6Þ overhead. MH pHðsjhpÞpppðhpÞ 2.BriefreviewofBayesruleandMCMC ThisprocedureisreferredtoastheMetropolis–Hastings(MH)algo- rithmaftertheoriginators[20,18].Thesameisdoneforeachofthe Givendatafromasystemofinterest,denotedbythevectorofN p=1,...,PparametersinhwhileholdingtheotherP(cid:3)1parameters observations s, and a model that describes those data w(h), the fixed.Thus,theapproachisreallydrawingsamplesfromthecondi- practitioner’s job is to estimate the model parameter vector tional posterior p ðh jh Þ. In other words, given that the other Hp p (cid:3)p h(cid:2)(h ,h ,...,h ).Typicallyonetakes‘‘good’’estimatestobethose parameters are fixed to their current values in the chain, sample 1 2 P that perform well in the face of the inevitable uncertainty (e.g. theposteriorforparameterp.ThisprocedureisreferredtoasGibbs noise)inthedata.Inthispaperitisassumedthattheuncertainty samplingandeliminatestheneedtoperformthehigh-dimensional takestheformofadditivenoisesuchthatourmeasureddataare integral required by Eq. (5). The above-described procedure Pleasecitethisarticleinpressas:NicholsJMetal.Bayesianidentificationofacrackedplateusingapopulation-basedMarkovChainMonteCarlomethod. ComputStruct(2011),doi:10.1016/j.compstruc.2011.03.013 J.M.Nicholsetal./ComputersandStructuresxxx(2011)xxx–xxx 3 constitutesoneiterationoftheMarkovchain.Startingwithinitial is the true posterior distribution. In what follows the authors use parameter values chosen from their respective prior distributions n=1todenotethetrueposteriordistributionfromwhichthede- hpð0Þ(cid:5)pppðhpÞthisprocess(accept/reject)repeatsKtimesyielding sired samples will be drawn. If an irreducible, a-periodic Markov thesamplesh (k),k=1...,K. chaincanbeconstructedthathasp (h(1:N))hasitsinvariantdistri- p C Thedistributionqðh(cid:4)jh Þcanbethoughtofasaruleforperturb- bution,thensamplesfromthemarginaldistribution p p inghp.Themagnitudeofthisperturbationwillbedecidedbythe Z parameters associated with q((cid:6)). For example, it is common to p ðhÞ¼ p ðhð1:NÞÞdhð2:NÞ H C choose the proposal distribution to be the Uniform distribution RN(cid:3)1 centeredathp,i.e. can be drawn where the notation RRN(cid:3)1dhð2:NÞ denotes the multi- qðh(cid:4)jh Þ¼Uðh (cid:3)A;h þAÞ ð7Þ dimensional integral over all parameters vectors other than p p p p h(cid:2)h(1). This is accomplished numerically by running ‘‘N’’ Markov where A is a real constant that plays the role of the perturbation chains concurrently, each exploring its own posterior distribution size. Typically A is tuned dynamically so that on average pn(h(n)). Each chain can be considered in turn, holding the others rMH=0.3(cid:3)0.5 (see [19] for a more thorough discussion and for fixed (Gibbs sampling), and samples from p1(h(1)) are retained as sample code). Heuristically it is easy to see how the approach samplesfromthedesiredjointposterior. works.Thealgorithmcontinuallyperturbseachparameter,checks Thereexistsagooddealoffreedominchoosingthepn((cid:6))andin whether or not a better fit (consistent with the prior) is achieved selectingdifferenttypesofproposalsforexploringtheparameter bycomputingr ,andkeepsparametervaluesthatdowellinthis space.Withregardtotheformeronewouldliketowiselychoose MH regard.Thefactthatintheendthechainofvaluesaresamplesfrom thepn((cid:6))soastofacilitateeasyexplorationoftheparameterspace. p ðh Þisactuallyquiteremarkable.However,onecanimmediately Whileotherapproachescanbeused(see[15]),inthispapertheso- Hp p seewhereproblemscanarise. calledtemperedsequenceofdistributionsisused Considerabimodalposteriordistributionofasingleparameter p ðhðnÞÞ(cid:2)pfnðhðnÞÞ ð10Þ h, n H pHðhÞ¼2pffi2ffi1ffipffiffiffiffifffiffirffiffiffi2ffiffie(cid:3)ðh(cid:3)lÞ2=ð2fr2Þ sfomraflnle2r(v0a,l1u]e.sFoorf ffn1=on1e, oofbctaoiunrsseu,cocnesesihvaeslythsemtorouteheprosvteerrsiioorn.sFoorf the original posterior. The idea is that the smoothed distributions 1 þ2pffi2ffiffipffiffiffiffiðffi1ffiffiffiffi(cid:3)ffiffiffiffiffifffiffiÞffiffirffiffiffi2ffiffie(cid:3)ðhþlÞ2=ð2ð1(cid:3)fÞr2Þ ð8Þ athreeyeaarseierstiflolrretlhaetierdctoorrtehsepotrnudeinpgosMtearrikoorvancdhatihnesretfooreexpstliolrlec,ayrreyt informationabouthighprobabilityregionsoftheparameterspace. forknown,positiveconstantsf,r,l.Theconstantf<1specifiesthe As an example consider N=4 separate instances of the bi-modal fractionofthedistributionvarianceassociatedwiththedistribution distribution(8)raisedtothef =1,0.5,0.1,0.0125powers,respec- peaksat+l,(cid:3)l.Forthisexamplethevaluesr=1,f =0.3,wereused n tively.ThesedistributionsareshowninFig.2.ThestandardMCMC alongwithtwodifferentvaluesforthepeaklocations,l=2,4.The algorithm will have an easier time exploring these smoother resulting distribution presents difficulties for the conventional parameter spaces. Of course only samples from the distribution MCMC algorithm for the l=4 case as shown in Fig. 1. For l=2 with f =1 are of interest, however the additional distributions thepeaksofthedistributionsarecloseenoughthatthealgorithm 1 canclearlyfacilitateefficientsamplingifonehasameansofpassing can easily move back and forth between the two regions of high information between chains. The chain associated with the true probability. Howeverfor l=4the MarkovChainquicklybecomes posteriorneedstobeinformedbythechainsexploringthesmooth- ‘‘trapped’’ in a single portion of the posterior distribution. This is erdistributions. simplyduetothefactthattheproposaldistributionisnotcapable Perhapsthesimplesttypeofmovetoaccomplishthisistheso- ofmovingthechaineasilyfromonepeaktotheother.Thesolution, called swap move. This is similar to a Metropolis–Hastings move however,isnotassimpleaschangingtheproposaldistributionto wheretheproposalistoconsiderswappingtheparametervalues giveuslargerperturbationstoourchain(i.e.increaseA).Ifthepro- in chains u, v. Assuming an equal probability of selecting chains posalvaluesareveryfarfromtheexistingvaluestheywillalmost u,vfromtheNpossibilities,thisswapisacceptedwithprobability alwaysberejected,thustheMarkovChainswilltakeanextremely longtimetoconverge.Anefficientsamplerisonethatwouldallow (p ðhðvÞÞp ðhðuÞÞ ) ustolocallyexploreahighprobabilityregionoftheposteriorwhile rswap¼min puðhðuÞÞpvðhðvÞÞ;1 : ð11Þ simultaneouslyprovideamechanismforcoveringlargedistancesin u v parameterspacetoreachotherhighprobabilityregions.Thepopu- Typicallyswapmovesarenotproposedaftereveryiterationinthe lation-basedMCMCapproach,describednext,wasdesignedspecif- MarkovChainsbutareperformedwithsomeprobability.Returning icallyforthisreason. tothebi-modalexample,considerN=4chainsexploringthecom- positetarget 3.Population-basedMCMC p ðhð1:4ÞÞ¼Y4 pfnðhðnÞÞ; ð12Þ C H AswithstandardMCMC,thegoalofthePop-MCMCalgorithmis n¼1 todrawsamplesfromsomedesiredposteriordistributionpH(h).A thusonehasfourMarkovchainsrunningconcurrentlywithf =1, very nice introduction to the topic of Population-based MCMC n 0.75,0.5,0.25,respectively.AftereachiterationintheMarkovchain (Pop-MCMC)isgivenbyJasraetal.[15].Thebasicideaistofirst aswapmoveisperformedwith50%probability.Thismoveconsists createanew,compositeposteriordensity ofuniformlyselectingtwoofthechainsandevaluatingEq.(11).If YN the proposal is accepted the values in the chains are exchanged p ðhð1:NÞÞ¼ p ðhðnÞÞ ð9Þ C n and the algorithm continues to the next iteration in the Markov n¼1 chains. Fig. 3 shows the results of this sampler All four Markov which is a function of the composite parameter vector chainsareinformingeachotherastothepresenceofmultiplepeaks h(1:N)=(h(1),...,h(N)). It is required that pn(h(n))=pH(h) for at least inthedistribution.Theendresultisthatthechainassociatedwith one n, i.e. one of the posteriors comprising the composite density f =1containssamplesfromthedesiredposteriordistribution;this 1 Pleasecitethisarticleinpressas:NicholsJMetal.Bayesianidentificationofacrackedplateusingapopulation-basedMarkovChainMonteCarlomethod. ComputStruct(2011),doi:10.1016/j.compstruc.2011.03.013 4 J.M.Nicholsetal./ComputersandStructuresxxx(2011)xxx–xxx (a) (b) 4 0.4 True PDF 3 Sampled PDF 2 0.3 1 Chain 1 −01 p(θ)θ 0.2 −2 0.1 −3 −4 −5 0 0 2000 4000 6000 8000 10000 −5 0 5 Iteration θ (c) (d) 6 1 True PDF Sampled PDF 5 0.8 4 0.6 1 n 3 θ) ai (θ h p C 0.4 2 0.2 1 0 0 0 2000 4000 6000 8000 10000 −5 0 5 Iteration θ Fig.1. MarkovChainsandestimatedbi-modalposteriordistributiongivenbyEq.(8)withl=2(a,b)andl=4(c,d),respectively. 1 canholdmultipleparametersfixedwhileallowingothergroupsof parameters to move simultaneously. Such a move is not allowed ζ=0.0125 withthestandardMetropolis-within-Gibbssampling(oneparame- 0.8 termovedatatime)orbytheswapmovewhichmovesallparame- ters simultaneously. A particularly effective move is to use 0.6 ζ=0.1 differentialevolutiontogenerateatrialvector.Differentialevolution θ) (DE)istheengineofapopularGAusedinsearchingcomplexparam- (n eterspace[21].Theapproachdrawsatrandomthreemembersofthe p 0.4 population,h(u),h(v),h(w)andgeneratesthetrialvector ζ=0.5 hðu0Þ¼hðuÞþcðhðvÞ(cid:3)hðwÞÞ ð13Þ 0.2 wherec isauser-definedconstant.EachofthePelementsinthis ζ=1 trial vector replace the elements of the original vector h(u) with 0 50% probability. This final step (keep new element or retain old) −10 −5 0 5 10 θ emulatesthe‘‘cross-over’’stepcommontomostGAs.Oncethetrial vectorhasbeengenerateditisaccepted/rejectedusingr (Eq.6). MH Fig.2. SuccessivepnðhÞ¼pfHnðhÞcorrespondingtothebi-modalposteriordistribu- However, in order to differentiate among the types of moves the tionEq.(8).Thetrueposteriordistributionfromwhichsamplesaresoughtisgiven acceptanceratiofortheDEmovewillbedenotedr . DE byf1=1. Itshouldmentionedthatothershaveproposedastraightcross- overmovewherebytheparametervectorissplitatapointm<P fortwo‘‘parent’’vectorsandgeneratingthetrialvectors: distributionisshowninFig.3bandcomparesfavorablytothetrue posteriordensity. hðu0Þ¼½hð1uÞ;...;hðmuÞ;h1ðvÞ;...;hPðvÞ(cid:7) ð14Þ Whiletheswapmoveisaneffectivemeansofcommunicating hðv0Þ¼½hðvÞ;...;hðvÞ;hðuÞ;...;hðuÞ(cid:7) 1 m 1 P among chains, a more sophisticated type of move is required for thestructuralestimationproblemconsiderednext.Particularlyfor and accepting with probability r (Eq. 11)[15]. This is indeed a swap multivariateparameterestimation,itisusefultodesignamovethat useful way of exchanging information between chains, however Pleasecitethisarticleinpressas:NicholsJMetal.Bayesianidentificationofacrackedplateusingapopulation-basedMarkovChainMonteCarlomethod. ComputStruct(2011),doi:10.1016/j.compstruc.2011.03.013 J.M.Nicholsetal./ComputersandStructuresxxx(2011)xxx–xxx 5 (a) 10 1)(k) 0 (θ −10 0 500 1000 1500 2000 10 (b) k) (2)θ( 0 0.5 True PDF −10 Sampled PDF 0 500 1000 1500 2000 0.4 10 (3)θ(k) 0 (θ)θ 0.3 −10 p 0 500 1000 1500 2000 0.2 10 0.1 k) 4)( 0 (θ −10 0 0 500 1000 1500 2000 −5 0 5 Iteration k θ Fig.3. (a)Markovchainsand(b)estimatedposteriordistributionassociatedwiththebi-modaldistributiongivenbyEq.(8). theauthorshavefoundthatthismechanismisalreadyprovidedfor isusedtodescribetheactualmotioninspaceandtime,w(x,y,t). inthedifferentialevolutionmove,thusthistypeofmoveisnotused The eigen-solution to this finite element model gives the modes inthispaper. intermsofthenodaldisplacements.Thevaluesattheinstrumen- The above-described moves clearly borrows from the genetic tationlocationsareinterpolatedviathenaturalneighbormethod. algorithm (GA) approach to optimization problems. In fact, Pop- Standardmodalanalysisproceedsfromthere.Thishybridsolution MCMCiseffectivelycombiningtheefficientsearchcapabilitiesof ismanytimesfasterthantimemarching,whileretainingtheflex- GAswiththepoweroftheBayesianMCMCapproachtosampling. ibilityinherentinthefiniteelementapproach. Bothtypesofmoves,swapanddifferentialevolution,willbeused Anin-housefiniteelementcodewasdevelopedforthispartic- in the structural dynamics example presented in subsequent ular application, so as to more easily wrap the MCMC process sections. (code)aroundthetimeseriessolution.Themodelparametersare thelocationofthecenterofthecrack(x ,y ),thecracklength crack crack (a),andtheorientationofthecrackmeasuredfromthepositivex- 4.Anefficientcrackedplatemodel axis(a).SeeFig.4c.Theelementsawayfromthecracktiparestan- dardeight-nodedquadrilateralMindlinserendipityelements,with Onecriticalcomponentofthisprocessisacomputationallyfast nodesateachcornerandatthemiddleofeachside.Adjacenttothe modelthatisalsoflexibleenoughtoincorporateavarietyofdam- cracktiptheeight-nodedquadsaremodifiedasdescribedinRefs. agemodes.Themodelpredictsthelateraldeflectionoftheplateat [23,24].Inshort,twocornernodesandthenodeinbetweenthem agivenpointinspaceandtime,w(x,y,t).Efficiencyoftheforward are moved to a single location, leaving a triangular element. The model is essential as the Pop-MCMC algorithm requires a very adjacent side nodes are then moved from the midpoints to one largenumberofmodel-to-datecomparisons,i.e.evaluater ,r , quarter of the side length. The result is an element with a stress MH swap rDE.Flexibilityisalsodesiredbecausewhilethecurrentstudycon- field that varies as p1ffirffi, where r is the distance from the collapsed siders damage consisting of only a single straight crack, it is the node(tobeplacedonthecracktip).SeeFig.4aandbforanillus- modelalonethatlimitsthetypeofdamagethatcanbeconsidered tration.ThisstressdistributionexactlymatchestheMitchellsolu- by this technique. A model that could handle multiple cracks, tion[25], which is the analytic stress field near a static crack tip. branched cracks, edge cracks, or corrosion would be valuable. A Because these augmented triangular elements capture the crack purely analytical model entailing virtually no computation time tip behavior, they are placed around the crack tip in a pinwheel wouldbeidealonbothcounts.Unfortunately,analyticalsolutions fashion. See Fig. 4d. This permits a more sparse mesh near the forthedynamicresponseofarectangularplatewithanarbitrary crack tip than would otherwise be necessary. This reduces the crack (position, orientation, and size) do not presently exist to numberofdegreesoffreedomrequiredand,hence,thecomputa- the knowledge of the authors. Solecki [22] produced a solution tiontime. for the naturalfrequencies ofa platewith asinglecrack, but did The cracked plate model involves several standard assump- not extend it to the myriad of other situations of future interest. tions. First, the model is linear and it is assumed that material A second possible family of models are finite element solutions. properties are known exactly (though these could be left as un- Theseclearlyhavethedesiredflexibility,butinvolvetimemarch- known parameters and found via the Bayesian/MCMC approach). ing and so are computationally intense. In this present work, a Thedeflectionsarepresumedsmallandcrackgrowthisnotcon- compromiseapproachistaken.Afiniteelementeigen-solutionis sidered;thelatterisreasonableundersmalldeflections.Thecrack usedtobuildinthefeaturesofthecracksingularity.Buildingoff is also assumed to remain open, such that impacts at the crack of these numerical frequencies and mode shapes, modal analysis interface are ignored. Any mass lost loss due to the crack is Pleasecitethisarticleinpressas:NicholsJMetal.Bayesianidentificationofacrackedplateusingapopulation-basedMarkovChainMonteCarlomethod. ComputStruct(2011),doi:10.1016/j.compstruc.2011.03.013 6 J.M.Nicholsetal./ComputersandStructuresxxx(2011)xxx–xxx Fig.4. (a)NodalconfigurationoftheMindlinserendipityelement,(b)resultinginterpolationfunctions,(c)schematicofplatewithcrackparametersand(d)resultingcourse FEMmesh. presumed negligible. Finally, there is only one crack and the tips Qðs;hÞ¼XT XM ðs (cid:3)w ðhÞÞ2; are at least one half of the crack length from the edges of the ij ij plate. Since finding smaller cracks is the more interesting and i j useful problem, the allowable zone is still a very large fraction thesequenceoftemperedposteriordistributionsisformedas oftheplatesarea.Eventakentogether,theauthorsdonotbelieve thaAtsthwesitehraenqyuifirenmiteenetlsemareenttooanoanlyesriosu,sc.onvergence of the mesh pnðhðnÞÞ¼pHðsjhðnÞÞfnppðhðnÞÞ¼ð2pprpð2hÞðfnnÞTÞM=2e(cid:3)2frn2nQðs;hðnÞÞ n must be verified. Due to the large number of iterations, it is not n¼1;...;N ð17Þ practicaltoverifyconvergenceforeveryperturbedsetofparame- ters.Instead,becauseallthemeshesusedarequalitativelysimilar, usingthesequencef1=1.0, convergence was checked on several representative parameter 1 vectors. fnþ1¼fn(cid:3)N; n¼1;...;N(cid:3)1 assuggestedbyJasraetal.[15].Theideahere,asinthetoyexam- 5.Implementation ple,istoexploreacompositeposteriorwherethe‘‘smoothed’’mar- ginal posteriors n=2,...,N are related to the true posterior of Inthispaperthestructureofinterestassumeaclampedrectan- interest,p (h(1)). gularplatemeasuring1.25m(cid:8)1mwiththicknessh=0.01mwith 1 This implementation of the pop-based algorithm proceeds as material properties E=209GPa (Young’s modulus), m=0.3 (Pois- follows.TheparametervectorsforeachoftheNchainsareinitial- son’sratio),andq=7850kg/m3(density).Itwillalsobeassumed izedbydrawingsamplesfromthepriors.Then,foreachiterationin thattheplatehasbeeninstrumentedwithMdisplacementsensors theMarkovchain,oneofthechainsn2Nisselectedwithuniform capableofsamplingtheplate’sresponsetoaninputat i=1,...,T probability,andastandardMHupdateisperformedforeachofthe equallyspacedpointsintime.Theobservedsignalmodelisthere- P parameters in h(n) using the Gibbs sampling strategy. Thus, for forewritten(asbefore) each parameter p=1,...,P one generates a candidate value using Eq.(7),evaluatestheratio s ¼w ðhÞþn ; i¼1;...T; j¼1;...;M ð15Þ ij ij ij r ¼Exp(cid:5)(cid:3) fn (cid:3)Qðs;h(cid:4)ðnÞjhðnÞÞ(cid:3)Qðs;hðnÞjhðnÞÞ(cid:4)(cid:6)ppðh(cid:4)pðnÞÞ rwahnedroemthpernoicjeasrse,tia.ek.eenacahsrnea(cid:5)lizNa(t0io,nrs2)o.fTahneimido,zdeerlo(-dmesecarnib,Gedauinsstiahne MH 2r2n p (cid:3)p (cid:3)p ppðhðpnÞÞ ij previoussection)isthereforeevaluatedatthesametimesi=1,...,T andacceptswithprobabilitymin(rMH,1).The(unknown)noisevar- and locations j=1,...,M as the observed data. Under this noise ianceassociatedwiththenthchain,rn,alsoneedstobesampled.It modelthelikelihoodfunctionforthedataisgivenby[19] hasbeendemonstrated[19]thatbychoosingavaguepriorforthis parameter,onemaydirectlysamplefromtheposteriorvia pHðsjhÞ¼ð2pr12ÞTM=2e(cid:3)2r12PTi PMj ðsij(cid:3)wijðhÞÞ2 ð16Þ r2n(cid:5)CðMT=2;2=1ðfnQðs;hðnÞÞÞ whereC(a,b)denotestheGammadistributionwithparametersa,b. Denotingthesum-squarederroroverallsensorsby Onceeachoftheparameters(includingthenoisevariance)havebeen Pleasecitethisarticleinpressas:NicholsJMetal.Bayesianidentificationofacrackedplateusingapopulation-basedMarkovChainMonteCarlomethod. ComputStruct(2011),doi:10.1016/j.compstruc.2011.03.013 J.M.Nicholsetal./ComputersandStructuresxxx(2011)xxx–xxx 7 updatedforchainn,twodifferentchainsu,v2[1,N]areselectedat thetrialvectorhðu0Þ.Theratioofpriorsintheacceptancecriteriais random with uniform probability. With 50% probability either a alsorequiredastheydonotcancelout.Onethereforeevaluates ‘‘swap’’moveor a‘‘DE’’ movebetweenthesechainsis performed. Fortheswapmove,thepriorratioswillcancel,thusoneevaluates r ¼Exp(cid:5)(cid:3) fu (cid:3)Qðs;hðu0ÞÞ(cid:3)Qðs;hðuÞÞ(cid:4)(cid:6)(cid:8)QPpppðhðpu0ÞÞ; ð19Þ r ¼Exp(cid:5)(cid:3)f (cid:7) 1 Qðs;hðvÞÞ(cid:3) 1 Qðs;hðuÞÞ(cid:8) DE 2r2u QPpppðhðpuÞÞ swap u 2r2 2r2 v u (cid:7) 1 1 (cid:8)(cid:6) andacceptthemovewithprobabilitymin(r ,1).Tosummarize,the (cid:3)fv 2r2Qðs;hðuÞÞ(cid:3)2r2Qðs;hðvÞÞ ; ð18Þ algorithmpicksachainatrandomandperfoDrEmsastandardMHup- u v date on each of the parameters, including the noise variance. The acceptingthemovewithprobabilitymin(r ,1).FortheDEmove algorithmthenselectschainsatrandomandperformseitheraswap swap onerequiresthreerandomlydrawnchains(seeEq.13)togenerate move (Eq. 18) or crossover (Eq. 19) with 50% probability. This 0.09 0.1 0.09 0.8 0.08 0.08 0.7 0.07 0.07 0.08 yc000...456 0000....00003456Q(s,x,y|a,θ)ccc a0.06 0000....00003456Q(s,a,θ|x,y)ccc 0.3 0.02 0.02 0.04 0.2 0.01 0.01 0.2 0.4 0.6 0.8 −50 0 50 x c θ c Fig.5. Mainargumentofthelikelihoodfunction,Q(s,h),plottedasafunctionof(a)cracklocationontheplateand(b)theparametersa(cracklength)andhc,crack orientation.Foreachplottheremainingparameterswereheldfixedattheirtruevalues.Trueparametervaluesaredenotedwithan‘‘X’’. (a) (b) 40 6 Prior Prior 35 Posterior Posterior 5 30 4 25 a) 20 x)c 3 p( p( 15 2 10 1 5 0 0 0 0.03 0.06 0.09 0.12 0.15 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 x a c (c) (d) 10 0.05 Prior Prior Posterior Posterior 8 0.04 6 0.03 y)c θ)c p( p( 4 0.02 2 0.01 0 0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 −90 −70 −50 −30 −10 10 30 50 70 90 y θ c c Fig.6. Estimatedposteriordensitiesassociatedwith(a)cracklength,(b,c)cracklocation,and(d)crackorientation.Thetrueparametervaluesarea=0.1,(xc,yc)=(0.3,0.5), andhc=30. Pleasecitethisarticleinpressas:NicholsJMetal.Bayesianidentificationofacrackedplateusingapopulation-basedMarkovChainMonteCarlomethod. ComputStruct(2011),doi:10.1016/j.compstruc.2011.03.013 8 J.M.Nicholsetal./ComputersandStructuresxxx(2011)xxx–xxx procedurerepeatsforsomenumberofiterationsuntilenoughsam- durationof2s(T=8000observations).Thefirstninemodeswere pleshavebeendrawnfromtheposteriordistribution. used in generating the solution, all well below the Nyquist fre- Thealgorithmisfairlysimpletoimplementinsoftware,how- quencyof2kHz.Thesignaltonoiseratiowassetat17dB(50:1) everisclearlycomputationallyintensive.Eachevaluationofaratio corresponding to the level of noise observed in previous (r ,r orr )requiressolvingtheforwardmodeldescribedin experiments. MH swap DE Section4.Theefficiencyoftheforwardmodelisthereforeofgreat Giventhisresponsethegoalistoestimatetheposteriordistri- importance in using MCMC in structural system identification butions associated with each of the crack parameters. Using the problems. standardMCMCalgorithmtheauthorswereunabletoconsistently identifytheparametersduetotheaforementionedproblemwith exploringcomplicatedlikelihoodfunctions.Fig.5showstheargu- 6.Crackedplateidentification mentofthelikelihood,Q(s,h)asafunctionofdifferentparameter combinations(holdingtheothersfixedattheirtruevalues)Aswith Therearefourparametersthatdeterminethestateofdamagein thesimple1-Dexample(Fig.1),onecanimmediatelyseethedif- theplatemodel:cracklengtha,cracklocation(x,y),andcrackori- ficulty presented by this estimation problem. Consider Fig. 5a. c c entation h. These damage parameters are fixed to the values The location parameter x has a minimum at the true location, c c a=0.1m, x =0.3m, y =0.5m, and h =30(cid:2). Further assume that x =0.3,however,thereisalsoa‘‘trough’’ofminimawithapartic- c c c c fourdisplacementsensorshavebeenplacedonthesurfaceofthe ularly low value at (x =0.85, y =0.28). A similarly complicated c c plateatthex–ylocations(0.375,0.375),(0.375,0.862),(0.862,0.37 likelihoodisshownfortheparameters a, h (Fig.5b).Theseplots c 5),(0.862,0.862).Theacquireddatawillconsistofeachsensor’sre- illustrate the fundamental difficulty in parameter identification sponse to four separate impacts (hammer strikes) at locations forstructuralsystems:multipledamagestatescanyieldverysimilar (0.29,0.275),(0.29,0.725),(0.96,0.275),and(0.96,0.275).Thesesen- structuralresponsedata. sorlocationswerechosensoastomaximizethesumofthefirstfour Forthisreasontheauthorshavegravitatedtowardthepopula- modesintheresponse.Noclaimsofoptimalityaremaderegarding tion-basedapproachtosamplingtheposterior.Intheexamplethat thischoice,infactfindingsensorlocationsthatproduceamorewell- follows,N=10chainswereusedforeachoftheparametersa,h,x, c c definedlikelihoodisanactiveareaofresearch. y.ThepriorforthecrackparameteraistaketobeaGammadis- c The data used in the identification procedure consist of four tribution,ppaðaÞ¼aa(cid:3)1ExCpð½a(cid:3)Þab=ab(cid:7)withparametersa=1.25,b=0.025. simulated impulse response signals, sampled at 4kHz for a Thus,essentiallynodamageisassumedattheoutset.Certainlythis (a) (b) 30 10 Prior Prior Posterior Posterior 25 8 20 6 a) 15 x)c p( p( 4 10 2 5 0 0 0 0.03 0.06 0.09 0.12 0.15 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 a x c (c) (d) 10 0.06 Prior Prior Posterior Posterior 0.05 8 0.04 6 y)c )θc0.03 p( p( 4 0.02 2 0.01 0 0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 −90 −70 −50 −30 −10 10 30 50 70 90 y θ c c Fig.7. Estimatedposteriordensitiesassociatedwith(a)cracklength,(b,c)cracklocation,and(d)crackorientation.Thetrueparametervaluesarea=0.1,(xc,yc)=(0.6,0.35), andhc=(cid:3)20. Pleasecitethisarticleinpressas:NicholsJMetal.Bayesianidentificationofacrackedplateusingapopulation-basedMarkovChainMonteCarlomethod. ComputStruct(2011),doi:10.1016/j.compstruc.2011.03.013 J.M.Nicholsetal./ComputersandStructuresxxx(2011)xxx–xxx 9 distribution can be altered to reflect a known level of damage, with different crack configurations. The authors view this as a however the more typical case is to assume the plate is healthy strength in difficult system identification problems, particularly andallowthedata(likelihood)todrivetheposterior.Theformof those in damage detection where quantifying uncertainty is crit- this prior (Gamma) was chosen based on the fact that the crack icaltomakinginformeddecisionsregardingthehealthofastruc- cannotbenegative.OnemightalsohavechosenaBetapriorwhich ture. Choosing to use more modes in the model helps the hasfinitesupportontheplate,howeverforthepriorparameters situation, however, in practice it becomes increasingly difficult chosen there is essentially no probability of the crack extending tomatchmodeltodataforhighermodes.Inthispaperthefocus offtheendoftheplate.Forthecracklocationparametersx,y uni- was on the first five modes of the response. c c formpriorswerechosenoverthespanofthecrack,i.e.x U(0,1.25), c y U(0,1) reflecting the fact that no a prior knowledge about the c 7.Conclusions cracklocationexists.Similarlytheprioroncrackanglewastaken ash U((cid:3)90,90)asingeneralonewillnotknowthecrackorienta- c This paper proposes a population-based Markov Chain Monte tionapriorieither.ThepriordistributionsaredisplayedinFig.6. Carlomethodforsolvingcomplicatedsystemidentificationprob- Initializing the Markov chains using these priors, the popula- lems in structural dynamics. Both the application and specific tion-basedMCMCalgorithmwasrun for60,000iterationswitha implementation of the method are new. The Pop-MCMC method burn-inof50,000iterations.Theremaining10,000iterationswere iswell-suitedtothetypesofposteriorparameterdistributionsof- storedandusedtoformtheposteriordensitiesshowninFig.6.In ten found in structural dynamics problems. These distributions this simulation, and in others the authors have looked at, crack tendtobemulti-modalandcaneasilyconfusestandardMCMCap- lengthisperhapsthemosteasilyidentifiedparameter.TheMarkov proaches. The population-based approach, by contrast, uses con- chainsinallpopulationstendtothetruevalue(a=0.1)afteronlya cepts from genetic algorithm search routines to more efficiently fewthousanditerations.Likewise,thelocationparametery exhib- c searchthisparameterspace.Thesamplingisdoneinsuchaway itsauni-modalposteriordistribution.Theotherparameters,x,h c c as to avoid becoming stuck in locally optimal solutions. The effi- are significantly more challenging to estimate as both show the cacy of the approach has been demonstrated in identifying crack clearpresenceofmultiplemaximainthelikelihood.Thelocation location,length,andorientationusingonlysimulatedimpulse-re- parameterx exhibitsmulti-modalbehavior.Inparticularonecan c sponse data from a cracked-plate model. Because the Pop-based seethemultipleminimaobservedinFig.5.Thetwopeakslocated MCMCmethodrequiresrepeatediterationsoftheforwardmodel, atx =0.3,x =0.4presentedasignificantchallengefortheregular c c great care was taken in developing an efficient model. This was MCMC algorithm. Because crack lengths were initially assumed accomplishedbyusingtailored‘‘serendipity’’elementstodescribe small, the first minimum encountered by the Markov chain was thestressfieldnearthecracktip.Thisallowsforgoodmodelcon- the one located at x =0.4 and the solution would often remain c vergencewithmanyfewerelementsthanareusedinstandardfi- here.Runningmultiplechainseasilyovercomesthisproblem.Per- nite element codes. The combination of efficient modeling and hapsthegreatestutilityoftheapproach,however,canbeseenin effective parameter identification routines can provide a wealth the identified orientation h. Regardless of crack geometry, there c ofinformationaboutthestate ofa structure.Boththeparameter arenearlyalwaysmultiplesolutionsforh thatcomeclosetomax- c estimates and the credible intervals associated with those esti- imizing the likelihood. The result are multiple well-defined local mates are obtained. This allows for confidence-based decisions maxima that easily trap the standard MCMC algorithm. By con- regardingthemaintenanceofastructureandalsoprovidesinfor- trast,Fig.6capturestherelativeheightsofthesemaximaindicat- mation needed in prognostics models for damage evolution in ingthatthehighestprobabilityforcrackorientationisatthetrue structures. value,h =30(cid:2). c Asasecondexample,Fig.7showstheidentifiedposteriordis- tributions for the case where the true crack parameters were set Acknowledgments to the values a=0.1, x =0.6, y =0.35, h =(cid:3)20. Additionally, the c c c Gammaprior(biasedtowardnocrack)waschangedtoaUniform TheauthorswouldliketoacknowledgetheOfficeofNavalRe- priorinordertodemonstratetheinsensitivityoftheapproachto search under Contract No. N00014-09-WX-2-1002 for providing thischoice.Itcouldbethatthepractitionerhasnoapriorinforma- fundingforthiswork. tionregardingthepresenceandlengthofacrack,thusauniform prior would be appropriate. All parameters are again correctly References identified, provided that the final estimate is taken as the maxi- mumaposteriorvalue.Again,oneseesmultiple‘‘good’’solutions [1] HoribeT,WatanabeK.Crackidentificationofplatesusinggeneticalgorithm. whichcanoftentrapthestandardMCMCalgorithm.Forexample, JSMEIntJ2006;49:403–10. the crack location parameter y has a fairly high probability of [2] StullCJ,EarlsCJ,AquinoW.Aposterioriinitialimperfectionidentificationin c shellbucklingproblems.ComputMethodsApplMechEng2008;198:260–8. being y =0.7 despite the fact that the true value is a factor of c [3] PanigrahiSK,ChakravertyS,MishraBK.Vibrationbaseddamagedetectionina twodifferent.Similarlymultiplepeaksfortheorientationparame- uniform strength beam using genetic algorithm. Meccanica 2009;44: terh canbeseen.TheauthorsnotethatthestandardMCMCalgo- 697–710. c [4] Hwang S-F, Wu J-C, He R-S. Identification of effective elastic constants of rithm was not able to sample from these complex posterior composite plates based on a hybrid genetic algorithm. Compos Struct distributions. 2009;90:217–24. The complexity of the posterior distributions is simply a con- [5] FonsecaJR,FriswellMI,MottersheadJE,LeesAW.Uncertaintyidentificationby sequence of trying to identify parameters that have little affect themaximumlikelihoodmethod.JSoundVib2005;288:587–99. [6] Xu GY, Zhu WD, Emory BH. Experimental and numerical investigation of on the global vibrations. It is obvious from the likelihood plots structuraldamagedetectionusingchangesinnaturalfrequencies.JVibAcoust (Fig. 5) that varying crack configurations can lead to nearly the 2007;129:686–700. samevibrationalresponse.Thisispartofthephysicsofstructural [7] KiureghianAD,DakessianT.Multipledesignpointsinfirstandsecond-order reliability.StructSafety1998;20:37–49. identification problems and the Bayesian approach, using the [8] Guo X, Bai W, Zhang W, Gao X. Confidence structural robust design and pop-MCMC algorithm, is correctly reflecting this uncertainty. optimization under stiffness and load uncertainties. Comput Methods Appl Rather than producing point estimates which may or may not MechEng2009;198:3378–99. [9] BeckJL,AuS-K.Bayesianupdatingofstructuralmodelsandreliabilityusing coincide with the true parameter value, the approach described Markov Chain Monte Carlo simulation. J Eng Mech – ASCE here allows for a full accounting of the probability associated 2002;128(4):380–91. Pleasecitethisarticleinpressas:NicholsJMetal.Bayesianidentificationofacrackedplateusingapopulation-basedMarkovChainMonteCarlomethod. ComputStruct(2011),doi:10.1016/j.compstruc.2011.03.013