Efficient Implementation Of Newton-Raphson Methods For Sequential Data Prediction

0.24 MB·

by  Burak C. Civek

1
Efficient Implementation Of Newton-Raphson
Methods For Sequential Data Prediction
Burak C. Civek and Suleyman S. Kozat, Senior Member, IEEE
Abstract—Weinvestigatetheproblemofsequentiallineardatapredictionforreallifebigdataapplications.Thesecondorder
algorithms,i.e.,Newton-RaphsonMethods,asymptoticallyachievetheperformanceofthe”best”possiblelineardatapredictormuch
fastercomparedtothefirstorderalgorithms,e.g.,OnlineGradientDescent.However,implementationofthesemethodsisnotusually
feasibleinbigdataapplicationsbecauseoftheextremelyhighcomputationalneeds.RegularimplementationoftheNewton-Raphson
7 MethodsrequiresacomputationalcomplexityintheorderofO(M2)foranM dimensionalfeaturevector,whilethefirstorder
1 algorithmsneedonlyO(M).Tothisend,inordertoeliminatethisgap,weintroduceahighlyefficientimplementationreducingthe
0 computationalcomplexityoftheNewton-RaphsonMethodsfromquadratictolinearscale.Thepresentedalgorithmprovidesthe
2 well-knownmeritsofthesecondordermethodswhileofferingthecomputationalcomplexityofO(M).Weutilizetheshiftednatureof
theconsecutivefeaturevectorsanddonotrelyonanystatisticalassumptions.Therefore,bothregularandfastimplementations
n
achievethesameperformanceinthesenseofmeansquareerror.Wedemonstratethecomputationalefficiencyofouralgorithmon
a
J reallifesequentialbigdatasets.Wealsoillustratethatthepresentedalgorithmisnumericallystable.

9
IndexTerms—Newton-Raphson,highlyefficient,bigdata,sequentialdataprediction.
1
✦

]
S
D
.
s 1 INTRODUCTION online learning algorithms in the data processingliterature
c
[ TECHNOLOGICAL developments in recent years have is the family of first order methods, e.g., Online Gradient
Descent [18], [19]. These methods only use the gradient
1 substantially increased the amount of data gathered
information to minimize the overall prediction cost. They
v from reallife systems[1], [2], [3], [4]. There exists a signifi-
8 cantdataflowthroughtherecentlyarisingapplicationssuch achieve logarithmic regret bounds that are theoretically
7 as large-scale sensor networks, information sensing mobile guaranteed to hold under certain assumptions [18]. Gradi-
3 entbasedmethodsarecomputationallymoreefficientcom-
devicesandwebbasedsocialnetworks[5],[6],[7].Thesize
5 pared to other families of online learning algorithms, i.e.,
as well as the dimensionality of this data strain the limits
1.0 mofacsusirvreenatmarocuhnittecotfurdeast.aSirnecseulptrionceasnsinegxcaensdsivsetorcionmgpsuutcah- wfohrearesexqtue∈ncReMof, tMhe-dciommepnusitoantiaolnfaelactuomrepvleexcittoyrsis{xont}lyt≥i0n,
0 the order of O(M). However, their convergence rates re-
tional cost, efficient machine learning and data processing
7
mainsignificantlyslowwhenachievinganoptimalsolution,
algorithmsareneeded[1],[8].
1
since no statistics other than the gradient is used [3], [16],
: Inthispaper,weinvestigatethewidelystudiedsequen-
v [19]. Even though gradient based methods suffer from this
tial prediction problem for high dimensional data streams.
i convergence issue, they are extensively used in big data
X Efficientpredictionalgorithmsspecifictobigdatasequences
applicationsduetosuchlowcomputationaldemand[20].
r havegreatimportanceforseveralreallifeapplicationssuch
a ashighfrequencytrading[9],forecasting[10],trendanalysis Different from the gradient based algorithms, the well
knownsecondorderNewton-Raphsonmethods,e.g,Online
[11], financial market [12] and locational tracking [13]. Un-
Newton Step, use the second order statistics, i.e., Hes-
fortunately,conventionalmethodsinmachinelearningand
sian of the cost function [18]. Hence, they asymptotically
dataprocessingliteraturesareinadequatetoefficientlyand
achieve the performance of the ”best” possible predictor
effectively process high dimensional data sequences [14],
much faster [17]. Existence of logarithmic regret bounds is
[15], [16]. Even though today’s computers have powerful
theoreticallyguaranteedforthisfamilyofalgorithmsaswell
processing units, traditional algorithms create a bottleneck
[18].Additionally,thesecondordermethodsarerobustand
evenforthatprocessingpowerwhenthedataisacquiredat
pronetohighlyvaryingdatastatistics,comparedtothefirst
highspeedsandtoolargeinsize[14],[15].
order methods, since they keep track of the second order
In order to mitigate the problem of excessive computa-
information[17],[21].Therefore,inthesenseofconvergence
tionalcost,weintroduce sequential,i.e.,online,processing,
rateandsteadystateerrorperformances,Newton-Raphson
where the data is neither stored nor reused, and avoid
methods considerably outperform the first order methods
”batch”processing.[16],[17].Onefamilyofthewellknown
[16], [17], [19]. However, the second order methods offer
a quadratic computational complexity, i.e., O(M2), while
• Burak C. Civek and S. S. Kozat are with the Department of Electrical
thegradientbasedalgorithmsprovidealinearrelation,i.e.,
andElectronicsEngineering,BilkentUniversity,Ankara06800, Turkey
(e-mail:{burak,kozat}@ee.bilkent.edu.tr). O(M). As a consequence, it is not usuallyfeasible for real-
life bigdata applicationsto utilizethe merits ofthe second
ThisworkissupportedinpartbyTurkishAcademyofSciencesOutstanding
ResearcherProgramme,TUBITAKContractNo.113E517. orderalgorithms[20].
SubmittedtoIEEETransactionsonKnowledgeandDataEngineering
In this paper, we study sequential data prediction, foranyn.Thedefinitionofwˆnisgivenfortheabsoluteloss
where the consecutive feature vectors are the shifted ver- case. To this end, we use the second order Online Newton
sions of each other, i.e., for a feature vector of xt = Step(ONS)algorithmtotraintheweightvectors.TheONS
[xt,xt 1,...,xt M]T, the upcomingvector is in the form of algorithm significantly outperforms the first order Online
xt+1 =− [xt+1,−xt,...,xt M+1]T. To this end, we introduce GradientDescent(OGD)algorithmintermsofconvergence
−
secondordermethodsforthisimportantproblemwithcom- rateandsteadystateerrorperformance sinceitkeepstrack
putationalcomplexityonlylinearinthedatadimension,i.e., ofthesecondorderstatisticsofthedatasequence[16],[18],
O(M).Weachievesuchanenormousreductionincomputa- [19].Theweightvectorateachtimeisupdatedas
tionalcomplexitysincethereareonlytwoentrieschanging
1
from xt to xt+1, where we avoid unnecessary calculations wt =wt−1− µA−t 1∇t, (2)
in each update. We do not use any statistical assumption
on the data sequence other than the shifted nature of the where µ isthe stepsizeand t M corresponds to
∈ R ∇ ∈ R
feature vectors. Therefore, we present an approach that is the gradient of the cost function ℓt(wt) at time t w.r.t. wt,
highlyappealingfor bigdata applicationssince it provides i.e., t , ℓt(wt). Here, the M M dimensional matrix
∇ ∇ ×
the merits of the Newton-Raphson methods with a much At isgivenby
lowercomputationalcost.
t
tialOdvaetraalpl,reinditchtiiosnpaaplgeor,riwthemintthraotdui)ceparnocoesnsleinseosnelqyuethne- At =Xi=0∇i∇Ti +αIM, (3)
currently available data without any storage, ii) efficiently
where α > 0 is chosen to guarantee that At is positive
implements the Newton-Raphsonmethods, i.e., the second
definite, i.e., At > 0, and hence, invertible. Selection of
ordermethodsiii)outperformsthegradientbasedmethods the parameters µ and α is crucial for good performance
intermsofperformance,iv)hasO(M)computationalcom-
[18]. Note that for the first order OGD algorithm, we have
plexity same as the first order methods and v) requires no
At = IM for all t, i.e., we do not use the second order
statistical assumptions on the data sequence. We illustrate
statisticsbutonlythegradientinformation.
the outstanding gains of our algorithm in terms of com- DefinitionofAtin(3)hasarecursivestructure,i.e.,At =
pduattaastieotsnaalndefficocmienpcayrebtyheurseisnugltitnwgoersreoqrupeenrtfioalrmreaanlcleisfewbitihg wAet−g1e+ta∇sttr∇aiTtg,hwtuitphdaanteinfriotimalAva−tlu11etoofAA−t−11u=sinαgItMhe.Hmeantrciex,
theregularNewton-Raphsonmethods. inversionlemma[22] −
A 1 TA 1
2 PROBLEM DESCRIPTION A−t 1 =A−t−11− 1−t+−1∇∇Ttt∇At−t−11∇−t−t1. (4)
Inthispaper,allvectorsarerealvaluedandcolumn-vectors. Multiplying both sides of (4) with t and inserting in (2)
Weuselowercase(uppercase)boldfaceletterstorepresent yields ∇
wvfoehrcettorherestvh(meecsatuotrrbicsxec.rsi)Tp.htTeihseiudoserenddtiinttyoarimyndatrtircaianxtseipstohrseaeptirtshedseeednnitometdeednbsyaiosInxMiTs, wt =wt−1− µ1(cid:20)1+A∇−tTt−1A1∇−t−t11∇t(cid:21). (5)
M M.WedenotetheM-dimensionalzerovectoras0M. Although the second order update algorithms provide
×
We study sequentialdata prediction, where we sequen- faster convergence rates and better steady state perfor-
tiallyobservearealvalueddatasequence xt t 0,xt . mances, computationalcomplexityissue prohibits theirus-
At each time t, after observing xt,xt 1,{...,x}t≥M+1 ∈, wRe ageinmostreallifeapplications[19],[22].Sinceeachupdate
glienneearramteoadneleastsimateofthedesire{ddata−,xˆt+1 ∈−R,us}inga imna(t4r)ixrewquitihresanthMemduilmtipelnicsaiotinoanlovfeacntoMr fo×rMxtdimenMsi,otnhael
xˆt+1 =wTt xt+ct, (1) computational complexity is in the order of O(∈M2R), while
thefirstorderalgorithmsjustneedO(M)multiplicationand
where xt M represents the feature vector of previous addition.Asanexample,inproteinstructureprediction,we
M sample∈s, iR.e., xt = [xt,xt 1,...,xt M+1]T. Here, wt have M = 1000 deeming the second order methods 1000
M and ct are the corr−espondin−g weight vector an∈d timesslowerthanthefirstorderOGDalgorithm[23].
R ∈ R
the offset variable respectively at time t. With an abuse of Inthenextsection,weintroduceasequentialprediction
notation, we combine the weight vector wt with the offset algorithm,which achievesthe performance ofthe Newton-
variablect, anddenote itbywt = [wt;ct], yieldingxˆt+1 = Raphson methods, while offering O(M) computational
wTt xt,wherext =[xt;1].Astheperformancecriterion,we complexitysameasthefirstordermethods.
use the widely studied instantaneous absolute loss as our
costfunction,i.e.,ℓt(wt)=ketk,wherethepredictionerror 3 EFFICIENT IMPLEMENTATION FOR COMPLEXITY
ateachtimeinstantisgivenbyet =xt+1−xˆt+1. REDUCTION
We adaptively learn the weight vector coefficients to
asymptoticallyachievethebestpossiblefixedweightvector Inthissection,weconstructanefficientimplementationthat
wˆn, which minimizes the total prediction error after n is based on the low rank property of the update matrices.
iteration,i.e., Insteadofdirectlyimplementingthesecondordermethods
asin(4)and(5),weuseunitaryandhyperbolictransforma-
n
wˆn =awrg∈mRMinXt=0kxt+1−wTxtk, tHioensssiaton-urepldaatetedtmheatwriexigAh−tt1v.ectorwt andtheinverseofthe
2
SubmittedtoIEEETransactionsonKnowledgeandDataEngineering
We work on time series data sequences, which directly of ηt. Moreover, we observe that the update term can be
implies that the feature vectors xt and xt+1 are highly expressedintermsofrank2matrices,whichisthekeypoint
related. More precisely, we have the following relation be- forthereductionofcomplexity.
tweenthesetwoconsecutivevectorsas Initially,weassumethatxt =0fort<0,whichdirectly
[xt+1,xTt ]=[xTt+1,xt−M+1]. (6) ifmoupnldiesasA−−11 = A−−12 = α1IM using (3). Therefore, ∆−1 is
1
Thisrelationshowsthatconsecutivedatavectorscarryquite ∆ 1 = diag 1,0,...,0, 1 . (16)
the same information, which is the basis of our algorithm. − α { − }
Weusetheinstantaneousabsoluteloss,whichisdefinedas Atthispoint,wedefinethe(M+1) 2dimensionalmatrix
×
ℓt(wt)=kxt+1−wTt xtk. (7) Λ−1 andthe2×2dimensionalmatrixΠ−1 as
1 1 0 ... 0 0 T 1 0
Adilctthioonugahptphleicaatbiosonlsu,tietliosssnoistwdiifdfeerlyenutisaebdleinwthheendeatta=pre0-. Λ−1 =rα(cid:20)0 0 ... 0 1(cid:21) ,Π−1 =(cid:20)0 −1(cid:21), (17)
However, we resolve this issue by setting a threshold ǫ toachievetheequalitygivenby
close to zeroand not updating the weight vector when the
absolute error is below this threshold, et < ǫ. From (4) ∆ 1 =Λ 1Π 1ΛT1. (18)
k k − − − −
and (5), the absolute loss results in the following update
rulesforwt andA−t 1, W2 persohpoewrt,yaitstahcehieenvdedo,fitthheodldisscfuosrsiaolnl,tthat0.oBncyeutshiengratnhke
wt =wt−1± µ1(cid:20)1+Ax−tTt−1A1x−tt11xt(cid:21), (8) greivfoernmiunl(a1t4io)nasofthedifference term,wer≥estatetheηt term
−
A−t 1 =A−t−11− A1−t+−11xxTttxATt−tA11−tx−1t1, (9) For the furthηetr=diηsct−us1s+ionx˜,TtwΛet−p1rΠetf−er1ΛmTt−at1rx˜ixt.notation a(1n9d)
−
since t = xt dependingonthesignoftheerror. represent(19)as
∇ ±
gtioprliIitcthaimtsioscnl,eeasAsre−ttnh1t1aiaxtltltyhaesrecsionumlt(sp8)le.frxoRitmaythotehfretthhmeaanstercgixoe-ntvtdiencgotordrmemarturaillx-- (cid:2)√ηt 0T2(cid:3)(cid:20)√0η2t(cid:21)=(cid:2)√ηt−1 x˜Tt Λt−1(cid:3)Θt−1(cid:20)Λ√Tt−ηt1−x˜1t(cid:21)(2,0)
A−t 11 from A−t−12 and then calculating the multiplication whereΘt 1 isdefinedas
aAAdnim−t−td−−−11e12cnxxosttmi−oi1pnnadawclivtimtihuda−oputurdaitlxalytaAena−tyrtu1ee1leax.,cphlwicihitteicrkhantiocowanl,cleuwdlageteedsoevfAetlh−to−e1p1Mxatd×firroeMmct We first −employ aΘtu−n1it,ary(cid:20)01G2ivΠe0ntT2s−1t(cid:21)ra.nsformation H(2G1,t)
Similarto[22],wefi−rstdefinethenormalizationtermof in order to zero out the second element of the vector
theupdaterulegivenin(8)as [√ηt 1,x˜Tt Λt 1] and then use a Θt 1-unitary Hyperbolic
ηt =1+xTt A−t 11xt. (10) rotat−ionHHB−,i.e.,HHB,tΘt−1HTHB−,t =Θt−1,toeliminate
− the last term [24]. Consequently, we achieve the following
Then, the difference between the consecutive terms ηt and updaterule
ηt 1 isgivenby
− ηt−ηt−1 =xTt A−t−11xt−xTt−1A−t−12xt−1. (11) whereHtre(cid:2)√prηetsen0tsT2t(cid:3)h=eo(cid:2)v√eηrat−ll1tranx˜sTtfoΛrtm−1a(cid:3)tiHontp,rocess.(E2x2-)
We define the (M + 1) 1 dimensional extended vector
xt˜hte=re[lxatti,oxnTt−gi1v]Tenainnd(6g)e,t×thefollowingtwoequalitiesusing fiTsohtreisnthcueepod(ftatt+eheg1si)evthetrsautnphsdefaontreemx,tawtnieoonramlmsaoalitnzraeicteiedosntihsteegrumuaprηdatan,ttheeodewdve[av2l2ue]re,.
ηt =1+x˜Tt (cid:20)A0−tT−11 00M(cid:21)x˜t, (12) ΛoftΛtetr−m1,iis.en.,otΛstu, fefixcpielincitt,lysi.nMceowreeovaelsr,oenveeendctahlecuinladtiinvgiduthael
M value of the vector A−t 11xt to update the weight vector
0 0T coefficients. −
ηt−1 =1+x˜Tt (cid:20)0M A−tM12(cid:21)x˜t. (13) We achieve the following equalities based on the same
− argumentthatweusedtoget(12)and(13)
Therefore,(11)becomes
ηt−ηt−1 =x˜Tt ∆t−1x˜t, (14) (cid:20)A−t−011xt(cid:21)=(cid:20)A0−tTM−11 00M(cid:21)x˜t, (23)
wheretheupdateterm∆t 1 isdefinedas
− 0 0 0T
∆t−1 ,(cid:20)A0−tTM−11 00M(cid:21)−(cid:20)00M A0−tTM12(cid:21). (15) (cid:20)A−t−12xt−1(cid:21)=(cid:20)0M A−tM−12(cid:21)x˜t. (24)
− Here,bysubtractingthesetwoequations,weget
Thisequationimpliesthatwe donot needthe exactvalues
tohfeAv−ta−l1u1eanodftAhe−t−d12efiinndeidviddiuffaelrleynacned∆itt−i1s fsourffithcieencatltcoulkantioown (cid:20)A−t−011xt(cid:21)=(cid:20)A−t−120xt−1(cid:21)+∆t−1x˜t. (25)
3
SubmittedtoIEEETransactionsonKnowledgeandDataEngineering
We emphasize that the same transformation Ht, which Algorithm 1:FastOnlineNewtonStep
aAwse−tf−o12ulxlsoetw−d1sttoo Age−tt−11√xηt,t,ifawlseoextrtaennsdfotrhmestraΛnts−f1ormtoedΛvtecatnodr 1 CDhaotao:s{exαt}1>t≥001,sweq0inudeo.n.wc.esiz0eM0 aTndthestepsizeµ;
Wwhee√dreeηn1two−t1ee√(cid:20)s(2Ahη6to−−t)w−11a20sxthBta−Ht1(cid:21)qt ==x˜BΛTt˜√Λt,1η−wtt1−[hx1eTtrAeHB−tt−1=1r,e0p(cid:20)]r√Teqsηaetnndts0QQtT2h(cid:21)e=,in(pΛ2u6t)t. 26345 ΛΠxw0−h=1i=x˜le=(cid:20)t010t=Mr≥[,x−αw00t1,(cid:20)d0(cid:21)x0o,=ΘTt ]T00=M; (cid:20),.01η.2.−10=0ΠT21(cid:21),1;ρ(cid:21)−1; =0M ;
matrixandB˜ statestheoutputmatrixofthetransformation. 7 xˆt+1 =wTt xt;
Then,thefollowingequalityisachieved 8 et =xt+1 xˆt+1;
−
sBinΘcet 1HBtT.iEsquΘattB˜i−n1Θgtt−uh1neBi˜etalTermy=,eBni.tesΘ.,ot−fB1mBHaTtrtiΘcet−s1inHbTtoBthTsi(d2e7=s) 109 BDet=erm(cid:20)√iρnηte0−t−a11(cid:21)Givx˜enΛTtsΛt−rto1−t1ati;onHG,t forB;
of(27−)yields 11 B´ =BHG,t;
12 DetermineaHyperbolicrotationHHB,t forB´;
q√ηt =(cid:20)A1−t−120xt−10(cid:21)+∆t−1x˜t,0 T 13 (cid:20)√ρ0ηt(cid:21)t 0ΛT2t=B´HHB,t;
qqT +QΠt−1QT = ηt 1 (cid:20)A−t 12xt 1(cid:21)(cid:20)A−t 12xt 1(cid:21) +∆t−1. 14 ifketk>ǫthen
− − − − − (28) 15 wt+1 =wt+sgn(et)µ1(cid:20)ρtη√ηt(cid:21);
t
We know from (25) that the left hand side of the first term
in(28)equalsto[xTt A−t−11,0]T andq isgivenby 1167 elsewt+1 =wt;
1 A 1 x 18 end
q = √ηt (cid:20) −t−01 t(cid:21). (29) 19 xt =[xt,xt−1,...,xt−M+1]T;
20 end
Hence, we identify the value of Qmatrix usingthe second
termin(28)as
0 0T Hence, we show that when the low rank property of the
M
QΠt−1QT =0M A−t−12xt−η1xTt−1A−t−12 dthifefeitreernacteiotnermt =∆it,isfoarchiieve0d. Tfohretre=foir−e,1t,hietitsrapnressfoerrvmeadtifoonr
 t 1  ≥
A 1 0 − 0 0 (30) in (26) gives all the necessary information and provides
+(cid:18)(cid:20) 0−tTM−1 0M(cid:21)−(cid:20)0TM A−tM12(cid:21)(cid:19) aupcdoamtepdleates update rule. As a result, the weight vector is
−
qqT,
wT(t1ohh5te)eh.rreeWefdoweirfkeef,neesroxeuwpnbacstnet−hidtaAuttht−tt−ihen12eg∆−tteth−rAim1s−tt−ee1ηq1rt1−umu1asAluiints−tyig−n1a2gtnhxidetts−uid1npxesdefiTta−rntt1iietnAirgoe−tn−lta1hg2teiieovvqneaun(la9uiln)es. whwetre=thewwintt−−d11iv,+iodstuhgaenlr(wvetai)slueµ1e(cid:20)Aof−t−ηA1t1−tx1t1(cid:21)x,tifiskefotkun>dǫby,mu(3lt3i-)
ofq yields plying (29) by √ηt, which is the l−eft upper most entry of
0 0 0 0 thetransformedmatrixB˜,andtakingthefirstM elements.
QΠt−1QT =(cid:18)(cid:20)0TM A−tM12(cid:21)−(cid:20)0TM A−tM11(cid:21)(cid:19) Tinhiteiacloizmaptiloentesaalngdorrietqhumiriesdpuropvdiadteeds.inAlgorithm1withall
+(cid:18)(cid:20)A0−tTM−11 00−M(cid:21)−(cid:20)00TM A0−tM−12(cid:21)(cid:19) 3, wThhiechprroescuesltssedinmthaetrcixomBphuatasttiohneadlicmoemnpsiloexnisty(Mof+O(2M)×).
A 1 0 A 1 0− Sincethereisnostatisticalassumptions,weobtainthesame
−t 1 M −t M
−(cid:18)(cid:20) 0T− 0 (cid:21)−(cid:20)0T 0 (cid:21)(cid:19) (31) errorratescomparedtotheregularimplementation.
M M
A 1 0 0 0
=(cid:20)0TM−t 0M(cid:21)−(cid:20)0TM A−tM11(cid:21) 4 SIMULATIONS
−
=∆ In this section, we illustrate the efficiency of our algorithm
t
=Λ Π ΛT. on real life sequential big datasets. We implement both
t t t the regular and the fast ONS algorithms on two different
ThisequalityimpliesthatΠistimeinvariant,i.e.,Πt 1 =Πt datasets,one of which is the CMU ARCTICspeechdataset
andQisgivenas − where a professional US English male speaker reads 1132
Q=Λt. (32) distinct utterances [25]. The recording is performed at 16
4
SubmittedtoIEEETransactionsonKnowledgeandDataEngineering
KHzandthereexistmorethan50millionsamples.Thesec-
onddatasetisaweatherforecastingtraceprovidedbyDavis
Weather station in Amherst, Massachusetts [26]. We also
implementthefirstorderOGDalgorithmonthisdatasetto
demonstrate the performance comparison with the second
order methods in terms of computational efficiency and
convergence rates. Temperature data was collected every
5 minutes during a period of 7 years from 2006 to 2013.
There exist more than 600 thousand samples. Hence, both
datasets are suitable for simulating big data scenarios. The
datasequencesarescaledtotherange[ 1,1].
−
4.1 ComputationalComplexityAnalysis
Asthefirstexperiment,weexaminethecomputationtimeof
both the proposedefficient ONS algorithm and the regular
ONS algorithm. We first work on the two partitions of
CMU ARCTIC speech dataset with lengths n = 5 107
and n = 2.5 107, and measure the corresponding·total
·
processing time. Sequences of different lengths are chosen
to illustrate the effect of increasing data length. For both
sequences, we choose feature vectors, xt M, with sev-
∈ R
eraldimensions rangingfrom M = 16toM = 128.InFig.
1.a, we demonstrate the computation time comparisons of
theregularandthefastimplementationsofONSalgorithm. Fig. 1: (a) Comparison of the computation time. R-ONS:
As expected from our results, complexity of the regularly Regular ONS, F-ONS: Fast ONS (b) Relative gain on the
implemented ONS algorithm shows a quadratic relation computationtimewithrespecttotheRegularONSandthe
with the dimension of the feature vectors, M. However, OGDalgorithmswhentheFastONSalgorithmisused.
the proposed efficient implementation provides a linear
relation.AsubstantialobservationfromFig.1.aisthat,with
anincreasingdimensionalityofthespaceoffeaturevectors, the regular and the fast ONS algorithms offer the same
the reduction in the complexity becomes outstanding. We error performances. However, there might occur negligible
alsoobservethatthegrowthinthedatasetlengthcausesthe numerical differences as a consequence of the finite pre-
samelineareffectonbothalgorithms,i.e.,doublingthetotal cision of real life computing systems. In the second part
lengthnresultsinthedoubledcomputationtime. of the experiments, we examine the effects of numerical
We then consider the weather forecasting temperature calculations on the mean square error curves of both the
dataset,whereinthiscasetotallengthofthedatasequence regular and the fast ONS algorithms. We first consider the
isnotasmuchasthepreviousdataset.Therefore,wespecif- CMUARCTICspeechdatasetwithn=5 104samplessince
·
ically concentrate on much largervalues for the dimension weobservethatthealgorithmsreachthesteadystateforthis
ofthefeaturevectors.Here,wechoosethedimensionofthe n.ThedimensionofthefeaturevectorsischosenasM =64,
spaceoffeaturevectorsrangingfromM =100toM =1000 and the learning rates are determined as 0.003 for both
and totallengthof the datasequence is n = 6 105. InFig. algorithms.Wedemonstratethecomparisonofmeansquare
·
1.b,weillustratetherelativegainoftheintroducedfastONS errorcurvesinFig.2.a.Adirectandsignificantobservation
algorithm with respect to the regular implementation of is that the efficient implementation is numerically stable.
ONSandtheOGDalgorithmintermsoftotalcomputation Thereisnoobservabledifference betweenthemeansquare
time.Weobservethattherelativecomputationtimegainof errorcurvesintermsofbothconvergenceandsteadystate.
the presented algorithm shows a significant improvement We work on the temperature tracking dataset as well
incomparisonwiththeregularONS,asthedatadimension forthenumericalstabilityanalysis.Inthiscase,we include
increases. However, we also observe that the relative gain the OGD algorithm into the comparison and increase the
falls intothe negativeregionwhen comparedwith thefirst dimensionoffeaturevectorstoM =400forallalgorithms.
order OGD algorithm. This is an expected result, since the Here, we examine the effect of large dimensionality on the
OGD uses only an M dimensional vector in each iteration, numericalperformance oftheproposedalgorithmandalso
whereas the fast ONS uses an (M + 2) 3 dimensional compare the error performances with the first order OGD
×
matrix and performs additional transformation operations algorithm. Similar to the first case, we only represent the
to update the weight vectors. Besides, the negative gain first 500 samples since the second order algorithms reach
remainsconstant,sincebothalgorithmseventuallyhavethe the steady state for this point. The learning rates are set
complexityinthesameorderofO(M). to 0.001 for the regular and the fast ONS and 0.1 for the
OGDalgorithm.InFig.2.b,weillustratethecorresponding
4.2 NumericalStabilityAnalysis
mean square error curves. Same as the previous analysis,
We theoretically show that the introduced algorithm effi- the fast ONS algorithm shows numerically no difference
ciently implements the ONS algorithm without any sta- compared to the regular ONS algorithm. Hence, even for
tistical assumptions or any information loss. Hence, both such high dimensional feature vectors, the proposed algo-
5
SubmittedtoIEEETransactionsonKnowledgeandDataEngineering
[2] C.Xu,Y.Zhang,R.Li,andX.Wu,“Onthefeasibilityofdistributed
kernelregressionforbigdata,”IEEETransactionsonKnowledgeand
DataEngineering,vol.28,no.11,pp.3041–3052,Nov2016.
[3] R.D’Ambrosio,W.Belhajali,andM.Barlaud,“Boostingstochastic
newtondescentforbigdatalargescaleclassification,”in2014IEEE
InternationalConferenceonBigData,Oct2014,pp.36–41.
[4] R.CouilletandM.Debbah,“Signalprocessinginlargesystems,”
IEEESignalProcessingMagazine,vol.24,pp.211–317,2013.
[5] R.Wolff,K.Bhaduri,andH.Kargupta,“Agenericlocalalgorithm
forminingdatastreamsinlargedistributedsystems,”IEEETrans-
actionsonKnowledgeandDataEngineering,vol.21,no.4,pp.465–
478,April2009.
[6] T. Wu, S. H. Yu, W. Liao, and C. S. Chang, “Temporal bipartite
projectionandlinkpredictionforonlinesocialnetworks,”in2014
IEEEInternationalConferenceonBigData,Oct2014,pp.52–59.
[7] Y. Yilmaz and X. Wang, “Sequential distributed detection in
energy-constrained wireless sensor networks,” IEEE Transactions
onSignalProcessing,vol.17,no.4,pp.335–339,2014.
[8] T.MoonandT.Weissman,“UniversalFIRMMSEfiltering,”IEEE
TransactionsonSignalProcessing,vol.57,no.3,pp.1068–1083,2009.
[9] R.Savani,“High-frequencytrading:Thefaster,thebetter?”IEEE
IntelligentSystems,vol.27,no.4,pp.70–73,July2012.
[10] P.GhoshandV.L.R.Chinthalapati,“Financialtimeseriesforecast-
ingusingagentbasedmodelsinequityandfxmarkets,”in2014
6thComputerScienceandElectronicEngineeringConference(CEEC),
Sept2014,pp.97–102.
[11] L. Deng, “Long-term trend in non-stationary time series with
nonlinear analysis techniques,” in 2013 6th InternationalCongress
Fig.2:(a)Datadimension:M = 64forbothalgorithms.(b)
onImageandSignalProcessing(CISP), vol.2,Dec2013, pp.1160–
Datadimension:M =400forallalgorithms. 1163.
[12] W. Cao, L. Cao, and Y. Song, “Coupled market behavior based
financialcrisisdetection,”inThe2013InternationalJointConference
onNeuralNetworks(IJCNN),Aug2013,pp.1–8.
rithmremainsnumericallystable.Additionally,weillustrate [13] Y. Ding, H. Tan, W. Luo, and L. M. Ni, “Exploring the use of
thatthefirstorderOGDalgorithmshowslessthanadequate diversereplicasforbiglocationtrackingdata,”in2014IEEE34th
performance intermsofconvergence rate.Therefore,nega- InternationalConferenceonDistributedComputingSystems(ICDCS),
June2014,pp.83–92.
tivegainonthecomputationtimeobservedintheprevious
[14] L. Bottou and Y. Le Cun, “On-line learning for very
experiment becomes insignificant when we consider the large data sets,” Applied Stochastic Models in Business and
meansquareerroranalysisinFig.2.b. Industry, vol. 21, no. 2, pp. 137–151, 2005. [Online]. Available:
http://dx.doi.org/10.1002/asmb.538
[15] L. Bottou and O. Bousquet, “The tradeoffs of large
scale learning,” in Advances in Neural Information
5 CONCLUSION Processing Systems, 2008, pp. 161–168. [Online]. Available:
http://leon.bottou.org/publications/pdf/nips-2007.pdf
In this paper, we investigate online sequentialdata predic- [16] N. Cesa-Bianchi and G. Lugosi, Prediction, Learning, and Games.
tion problem for high dimensional data sequences. Even Cambridge:CambridgeUniversityPress,2006.
[17] A. C. Singer, S. S. Kozat, and M. Feder, “Universal linear least
thoughthesecondorderNewton-Raphsonmethodsachieve
squares prediction: upper and lower bounds,” IEEE Transactions
superiorperformance,comparedtothegradientbasedalgo-
onInformationTheory,vol.48,no.8,pp.2354–2362,Aug2002.
rithms, the problem of extremely high computational cost [18] E.Hazan,A.Agarwal,andS.Kale,“Logarithmicregretalgorithms
prohibits their usage in real life big data applications. For foronlineconvexoptimization,”MachineLearning,vol.69,no.2-3,
pp.169–192,2007.
an M dimensional feature vector, the computational com-
[19] D.BertsimasandJ.N.Tsitsiklis,Introductiontolinearoptimization,
plexity of these methods increases in the order of O(M2). ser. Athena scientific series in optimization and neural
Tothisend,weintroduceahighlyefficientimplementation computation. Belmont(Mass.):AthenaScientific,1997.[Online].
that reduces the computational complexity of the Newton- Available:http://opac.inria.fr/record=b1094316
Raphson methods from O(M2) to O(M). The presented [20] E. K. P. Chong and S. H. Zak, An introduction to optimization,
ser. Wiley-Interscience series in discrete mathematics and
algorithm does not require any statistical assumption on optimization. New York: Wiley, 2008. [Online]. Available:
the data sequence. We only use the similarity between the http://opac.inria.fr/record=b1128546
[21] S.S.Kozat,A.T.Erdogan,A.C.Singer,andA.H.Sayed,“Steady-
consecutive feature vectors without any information loss.
statemseperformanceanalysisofmixtureapproachestoadaptive
Hence,ouralgorithmofferstheoutstandingperformanceof filtering,”IEEETransactionsonSignalProcessing,vol.58,no.8,pp.
the secondorder methods withthe low computationalcost 4050–4063,Aug2010.
of the first order methods. We illustrate that the efficient [22] A.H.Sayed,FundamentalsofAdaptiveFiltering. NJ:JohnWiley&
Sons,2003.
implementation of Newton-Raphson methods attains sig-
[23] J.Cheng, A.N. Tegge,and P. Baldi, “Machine learning methods
nificant computationalgains,asthe datadimension grows. forproteinstructureprediction,”IEEEReviewsinBiomedicalEngi-
Wealsoshowthatouralgorithmisnumericallystable. neering,vol.1,pp.41–49,2008.
[24] A.H.Sayed,AdaptiveFilters. NJ:JohnWiley&Sons,2008.
[25] J. Kominek and A. W. Black, “Cmu arctic databases.” [Online].
Available:http://www.festvox.org/cmu arctic/index.html
REFERENCES [26] M. Liberatore, “Umass trace repository.” [Online]. Available:
http://traces.cs.umass.edu/index.php/Sensors/Sensors
[1] X. Wu, X. Zhu, G. Q. Wu, and W. Ding, “Data mining with big
data,”IEEETransactionsonKnowledgeandDataEngineering,vol.26,
no.1,pp.97–107,Jan2014.
6