首页 Python的隐马尔科夫HMMLearn库的应用教学(1)

Python的隐马尔科夫HMMLearn库的应用教学(1)

举报
开通vip

Python的隐马尔科夫HMMLearn库的应用教学(1)PythonHMMLearnTutorialEditedBy毛片物语hmmlearn implementstheHiddenMarkovModels(HMMs).TheHMMisagenerativeprobabilisticmodel,inwhichasequenceofobservable \(\mathbf{X}\) variablesisgeneratedbyasequenceofinternalhiddenstates \(\mathbf{Z}\).Thehiddenstatesarenotbeobser...

Python的隐马尔科夫HMMLearn库的应用教学(1)
PythonHMMLearnTutorialEditedBy毛片物语hmmlearn implementstheHiddenMarkovModels(HMMs).TheHMMisagenerativeprobabilisticmodel,inwhichasequenceofobservable \(\mathbf{X}\) variablesisgeneratedbyasequenceofinternalhiddenstates \(\mathbf{Z}\).Thehiddenstatesarenotbeobserveddirectly.Thetransitionsbetweenhiddenstatesareassumedtohavetheformofa(first-order)Markovchain.Theycanbespecifiedbythestartprobabilityvector \(\boldsymbol{\pi}\) andatransitionprobabilitymatrix \(\mathbf{A}\).Theemissionprobabilityofanobservablecanbeanydistributionwithparameters \(\boldsymbol{\theta}\) conditionedonthecurrenthiddenstate.TheHMMiscompletelydeterminedby \(\boldsymbol{\pi}\), \(\mathbf{A}\) and \(\boldsymbol{\theta}\).TherearethreefundamentalproblemsforHMMs:●Giventhemodelparametersandobserveddata,estimatetheoptimalsequenceofhiddenstates.●Giventhemodelparametersandobserveddata,calculatethelikelihoodofthedata.●Givenjusttheobserveddata,estimatethemodelparameters.ThefirstandthesecondproblemcanbesolvedbythedynamicprogrammingalgorithmsknownastheViterbialgorithmandtheForward-Backwardalgorithm,respectively.ThelastonecanbesolvedbyaniterativeExpectation-Maximization(EM)algorithm,knownastheBaum-Welchalgorithm.References:[Rabiner89]LawrenceR.Rabiner“AtutorialonhiddenMarkovmodelsandselectedapplicationsinspeechrecognition”,ProceedingsoftheIEEE77.2,pp.257-286,1989.[Bilmes98]JeffA.Bilmes,“AgentletutorialoftheEMalgorithmanditsapplicationtoparameterestimationforGaussianmixtureandhiddenMarkovmodels.”,1998.Availablemodelshmm.GaussianHMMHiddenMarkovModelwithGaussianemissions.hmm.GMMHMMHiddenMarkovModelwithGaussianmixtureemissions.hmm.MultinomialHMMHiddenMarkovModelwithmultinomial(discrete)emissionsReadon fordetailsonhowtoimplementanHMMwithacustomemissionprobability.BuildingHMMandgeneratingsamplesYoucanbuildanHMMinstancebypassingtheparametersdescribedabovetotheconstructor.Then,youcangeneratesamplesfromtheHMMbycalling sample.>>>importnumpyasnp>>>fromhmmlearnimporthmm>>>np.random.seed(42)>>>model=hmm.GaussianHMM(n_components=3,covariance_type="full")>>>model.startprob_=np.array([0.6,0.3,0.1])>>>model.transmat_=np.array([[0.7,0.2,0.1],...              [0.3,0.5,0.2],...              [0.3,0.3,0.4]])>>>model.means_=np.array([[0.0,0.0],[3.0,-3.0],[5.0,10.0]])>>>model.covars_=np.tile(np.identity(2),(3,1,1))>>>X,Z=model.sample(100)Thetransitionprobabilitymatrixneednottobeergodic.Forinstance,aleft-rightHMMcanbedefinedasfollows:>>>lr=hmm.GaussianHMM(n_components=3,covariance_type="diag",...          init_params="cm",params="cmt")>>>lr.startprob_=np.array([1.0,0.0,0.0])>>>lr.transmat_=np.array([[0.5,0.5,0.0],...            [0.0,0.5,0.5],...            [0.0,0.0,1.0]])Ifanyoftherequiredparametersaremissing, sample willraiseanexception:>>>hmm.GaussianHMM(n_components=3)>>>X,Z=model.sample(100)Traceback(mostrecentcalllast):...sklearn.utils.validation.NotFittedError:ThisGaussianHMMinstanceisnotfittedyet.Call'fit'withappropriateargumentsbeforeusingthismethod.FixingparametersEachHMMparameterhasacharactercodewhichcanbeusedtocustomizeitsinitializationandestimation.EMalgorithmneedsastartingpointtoproceed,thuspriortotrainingeachparameterisassignedavalueeitherrandomorcomputedfromthedata.Itispossibletohookintothisprocessandprovideastartingpointexplicitly.Todoso1.ensurethatthecharactercodefortheparameterismissingfrom init_params andthen2.settheparametertothedesiredvalue.Forexample,consideranHMMwithexplicitlyinitializedtransitionprobabilitymatrix>>>model=hmm.GaussianHMM(n_components=3,n_iter=100,init_params="mcs")>>>model.transmat_=np.array([[0.7,0.2,0.1],...              [0.3,0.5,0.2],...              [0.3,0.3,0.4]])Asimilartrickappliestoparameterestimation.Ifyouwanttofixsomeparameterataspecificvalue,removethecorrespondingcharacterfrom params andsettheparametervaluebeforetraining.Examples:●SamplingfromHMMTrainingHMMparametersandinferringthehiddenstatesYoucantrainanHMMbycallingthe fit method.Theinputisamatrixofconcatenatedsequencesofobservations(aka samples)alongwiththelengthsofthesequences(see Workingwithmultiplesequences).Note,sincetheEMalgorithmisagradient-basedoptimizationmethod,itwillgenerallygetstuckinlocaloptima.Youshouldingeneraltrytorun fit withvariousinitializationsandselectthehighestscoredmodel.Thescoreofthemodelcanbecalculatedbythe score method.Theinferredoptimalhiddenstatescanbeobtainedbycalling predict method.The predict methodcanbespecifiedwithdecoderalgorithm.CurrentlytheViterbialgorithm("viterbi"),andmaximumaposterioriestimation("map")aresupported.Thistime,theinputisasinglesequenceofobservedvalues.Note,thestatesin remodel willhaveadifferentorderthanthoseinthegeneratingmodel.>>>remodel=hmm.GaussianHMM(n_components=3,covariance_type="full",n_iter=100)>>>remodel.fit(X) GaussianHMM(algorithm='viterbi',...>>>Z2=remodel.predict(X)MonitoringconvergenceThenumberofEMalgorithmiterationisupperboundedbythe n_iter parameter.Thetrainingproceedsuntil n_iter stepswereperformedorthechangeinscoreislowerthanthespecifiedthreshold tol.Note,thatdependingonthedataEMalgorithmmayormaynotachieveconvergenceinthegivennumberofsteps.Youcanusethe monitor_ attributetodiagnoseconvergence:>>>remodel.monitor_ ConvergenceMonitor(history=[...],iter=12,n_iter=100,tol=0.01,verbose=False)>>>remodel.monitor_.convergedTrueWorkingwithmultiplesequencesAlloftheexamplessofarwereusingasinglesequenceofobservations.Theinputformatinthecaseofmultiplesequencesisabitinvolvedandisbestunderstoodbyexample.Considertwo1Dsequences:>>>X1=[[0.5],[1.0],[-1.0],[0.42],[0.24]]>>>X2=[[2.4],[4.2],[0.5],[-0.24]]Topassbothsequencesto fit or predict firstconcatenatethemintoasinglearrayandthencomputeanarrayofsequencelengths:>>>X=np.concatenate([X1,X2])>>>lengths=[len(X1),len(X2)]Finallyjustcallthedesiredmethodwith X and lengths:>>>hmm.GaussianHMM(n_components=3).fit(X,lengths) GaussianHMM(algorithm='viterbi',...Examples:●GaussianHMMofstockdataSavingandloadingHMMAftertraninganHMMcanbeeasilypersistedforfutureusewiththestandard pickle moduleoritsmoreefficientreplacementinthe joblib package:>>>fromsklearn.externalsimportjoblib>>>joblib.dump(remodel,"filename.pkl")["filename.pkl"]>>>joblib.load("filename.pkl") GaussianHMM(algorithm='viterbi',...ImplementingHMMswithcustomemissionprobabilitiesIfyouwanttoimplementotheremissionprobability(e.g.Poisson),youhavetosubclass _BaseHMMandoverridethefollowingmethodsbase._BaseHMM._init(X,lengths)Initializesmodelparameterspriortofitting.base._BaseHMM._check()Validatesmodelparameterspriortofitting.base._BaseHMM._generate_sample_from_state(state)Generatesarandomsamplefromagivencomponent.base._BaseHMM._compute_log_likelihood(X)Computesper-componentlogprobabilityunderthemodel.base._BaseHMM._initialize_sufficient_statistics()InitializessufficientstatisticsrequiredforM-step.base._BaseHMM._accumulate_sufficient_statistics(...)Updatessufficientstatisticsfromagivensample.base._BaseHMM._do_mstep(stats)PerformstheM-stepofEMalgorithm.SamplingfromHMMThisscriptshowshowtosamplepointsfromaHidenMarkovModel(HMM):weusea4-componentswithspecifiedmeanandcovariance.Theplotshowthesequenceofobservationsgeneratedwiththetransitionsbetweenthem.Wecanseethat,asspecifiedbyourtransitionmatrix,therearenotransitionbetweencomponent1and3.print(__doc__)importnumpyasnpimportmatplotlib.pyplotaspltfromhmmlearnimporthmmPrepareparametersfora4-componentsHMMInitialpopulationprobabilitystartprob=np.array([0.6,0.3,0.1,0.0])#Thetransitionmatrix,notethattherearenotransitionspossible#betweencomponent1and3transmat=np.array([[0.7,0.2,0.0,0.1],[0.3,0.5,0.2,0.0],[0.0,0.3,0.5,0.2],[0.2,0.0,0.2,0.6]])#Themeansofeachcomponentmeans=np.array([[0.0, 0.0],[0.0,11.0],[9.0,10.0],[11.0,-1.0]])#Thecovarianceofeachcomponentcovars=.5*np.tile(np.identity(2),(4,1,1))#BuildanHMMinstanceandsetparametersmodel=hmm.GaussianHMM(n_components=4,covariance_type="full")#Insteadoffittingitfromthedata,wedirectlysettheestimated#parameters,themeansandcovarianceofthecomponentsmodel.startprob_=startprobmodel.transmat_=transmatmodel.means_=meansmodel.covars_=covars#GeneratesamplesX,Z=model.sample(500)#Plotthesampleddataplt.plot(X[:,0],X[:,1],".-",label="observations",ms=6,mfc="orange",alpha=0.7)#Indicatethecomponentnumbersfori,minenumerate(means):plt.text(m[0],m[1],'Component%i'%(i1),size=17,horizontalalignment='center',bbox=dict(alpha=.7,facecolor='w'))plt.legend(loc='best')plt.show()Totalrunningtimeofthescript: (0minutes0.528seconds)GaussianHMMofstockdataThisscriptshowshowtouseGaussianHMMonstockpricedatafromYahoo!finance.Formoreinformationonhowtovisualizestockpriceswithmatplotlib,pleasereferto date_demo1.py ofmatplotlib.from__future__importprint_functionimportdatetimeimportnumpyasnpfrommatplotlibimportcm,pyplotaspltfrommatplotlib.datesimportYearLocator,MonthLocatortry:frommatplotlib.financeimportquotes_historical_yahoo_ochlexceptImportError:#ForMatplotlibpriorto1.5.frommatplotlib.financeimport(quotes_historical_yahooasquotes_historical_yahoo_ochl)fromhmmlearn.hmmimportGaussianHMMprint(__doc__)GetquotesfromYahoo!financequotes=quotes_historical_yahoo_ochl("INTC",datetime.date(1995,1,1),datetime.date(2012,1,6))#Unpackquotesdates=np.array([q[0]forqinquotes],dtype=int)close_v=np.array([q[2]forqinquotes])volume=np.array([q[5]forqinquotes])[1:]#Takediffofclosevalue.Notethatthismakes#``len(diff)=len(close_t)-1``,therefore,otherquantitiesalso#needtobeshiftedby1.diff=np.diff(close_v)dates=dates[1:]close_v=close_v[1:]#Packdiffandvolumefortraining.X=np.column_stack([diff,volume])RunGaussianHMMprint("fittingtoHMManddecoding...",end="")#MakeanHMMinstanceandexecutefitmodel=GaussianHMM(n_components=4,covariance_type="diag",n_iter=1000).fit(X)#Predicttheoptimalsequenceofinternalhiddenstatehidden_states=model.predict(X)print("done")Out:fittingtoHMManddecoding...donePrinttrainedparametersandplotprint("Transitionmatrix")print(model.transmat_)print()print("Meansandvarsofeachhiddenstate")foriinrange(model.n_components):print("{0}thhiddenstate".format(i))print("mean=",model.means_[i])print("var=",np.diag(model.covars_[i]))print()fig,axs=plt.subplots(model.n_components,sharex=True,sharey=True)colours=cm.rainbow(np.linspace(0,1,model.n_components))fori,(ax,colour)inenumerate(zip(axs,colours)):#Usefancyindexingtoplotdataineachstate.mask=hidden_states==iax.plot_date(dates[mask],close_v[mask],".-",c=colour)ax.set_title("{0}thhiddenstate".format(i))#Formattheticks.ax.xaxis.set_major_locator(YearLocator())ax.xaxis.set_minor_locator(MonthLocator())ax.grid(True)plt.show()Out:Transitionmatrix[[ 9.79217702e-01 3.55338063e-15 2.72110180e-03 1.80611963e-02][ 1.21602143e-12 7.73505488e-01 1.85141936e-01 4.13525763e-02][ 3.25253603e-03 1.12652335e-01 8.83404334e-01 6.90794633e-04][ 1.18928464e-01 4.20116465e-01 1.91329669e-18 4.60955072e-01]]Meansandvarsofeachhiddenstate0thhiddenstatemean= [ 2.40689227e-02 4.97390967e07]var= [ 7.42026137e-01 2.49469027e14]1thhiddenstatemean= [ 2.19283454e-02 8.82098779e07]var= [ 1.26266869e-01 5.64899722e14]2thhiddenstatemean= [ 7.93313395e-03 5.43199848e07]var= [ 5.34313422e-02 1.54645172e14]3thhiddenstatemean= [-3.64907452e-01 1.53097324e08]var= [ 2.72118688e00 5.88892979e15]Totalrunningtimeofthescript: (0minutes2.205seconds)
本文档为【Python的隐马尔科夫HMMLearn库的应用教学(1)】,请使用软件OFFICE或WPS软件打开。作品中的文字与图均可以修改和编辑, 图片更改请在作品中右键图片并更换,文字修改请直接点击文字进行修改,也可以新增和删除文档中的内容。
该文档来自用户分享,如有侵权行为请发邮件ishare@vip.sina.com联系网站客服,我们会及时删除。
[版权声明] 本站所有资料为用户分享产生,若发现您的权利被侵害,请联系客服邮件isharekefu@iask.cn,我们尽快处理。
本作品所展示的图片、画像、字体、音乐的版权可能需版权方额外授权,请谨慎使用。
网站提供的党政主题相关内容(国旗、国徽、党徽..)目的在于配合国家政策宣传,仅限个人学习分享使用,禁止用于任何广告和商用目的。
下载需要: 免费 已有0 人下载
最新资料
资料动态
专题动态
is_751406
暂无简介~
格式:doc
大小:57KB
软件:Word
页数:22
分类:
上传时间:2022-08-01
浏览量:1