molecular 
 dynamics ( md ) has had a long and important history in 
 computational biology and has made many recent contributions 
 toward understanding a variety of biological phenomena , including 
 protein folding , protein ligand binding , and conformational change 
 of large biomolecules .
though running a large - scale simulation has historically been a 
 demanding task , recent advancements in both software and hardware 
 have made simulations of complex systems possible .
nonetheless , because of the high - dimensional nature of the resulting 
 time series , the generation of many md trajectories is not sufficient 
 for providing scientific insight .
in fact , data analysis has increasingly 
 become a major barrier in the application of md simulations .
one analysis technique attempts to estimate 
 the slowly equilibrating 
 collective motions , or eigenprocesses , contained in a simulation ( e.g. , 
 proteins folding from unfolded conformations , ligands binding to a 
 target , enzymes shifting from an inactive to active state , etc . ) . 

these eigenprocesses are useful as they can be measured experimentally , 
 but by constructing them from a simulation , we can provide a level 
 of detail that is unavailable to experiment alone .
markov state models 
 ( msms ) are one such method that has been useful in a variety of contexts .
briefly , an msm describes 
 continuous dynamics in phase space as a jump process over a finite 
 number of states that partition the space .
practically , the construction of 
 an msm is an active area of research 
 as there are significant hurdles that must be overcome in order to 
 build a model that accurately describes the dynamics contained in 
 a simulation .
for example , an msm requires 
 the definition of a state decomposition  that is , a partition 
 of phase space into discrete states .
the accuracy of the model relies 
 heavily on the quality of this state decomposition , but although there 
 have been a number of improvements to this process , there is no consensus on the best method for defining a state space .
an msm , however , is merely one way to 
 estimate a system s 
 eigenprocesses from a simulation .
therefore , we can hope to overcome 
 the practical challenges faced when constructing an msm by approximating 
 the dynamics within a different scheme . toward this end
, we introduce 
 an extension of time - structure based independent components analysis 
 ( tica ) using the kernel trick for estimating the characteristic eigenprocesses 
 contained in a simulation without using an msm . first , 
 we lay the theoretical foundation for the method , then formulate the 
 problem and derive its solution , and finally , we use the method to 
 analyze several toy models as well as md simulations of the fip35 
 ww domain from shaw et al .
in this 
 section , we briefly discuss markovian dynamics and the transfer operator 
 approach but leave most details to other works .
we also note that , for clarity , most of our notation is borrowed 
 from prinz et al .
we are interested 
 in discrete - time dynamics in some phase space  that are markovian , 
 ergodic , and time - reversible .
this density is properly normalized such 
 that  dx (x ) = 1 . let x(t ) denote a markov process 
 in . because of the markov property , there exists a transition 
 probability density   such that p(x , y;)dy is the probability 
 of being in some volume ( dy ) around y after time  given that the system started at x. the dynamical progression of the markov chain can 
 be described 
 via a propagator operator , which maps a probability density at time t , pt(x ) , to a new density at time t +  
 later , pt+(x):1additionally , there exists a corresponding 
 operator called the transfer operator defined such that2 the functions ut(x ) are related to probability densities by3(while the 
 propagator provides a more intuitive 
 operation , the transfer operator s eigenfunctions will be the 
 eigenprocesses we actually estimate in the following methods . )
the action of the transfer operator can be decomposed into a series 
 of relaxation processes , which correspond to its eigenfunctions , i(x ) , and their eigenvalues , 
 i:4where5 the eigenfunctions of the 
 transfer operator are orthogonal and 
 normalized such that6as a result 
 of this decomposition , we can 
 model the dynamics of the markov chain directly if we have an estimate 
 of the eigenfunctions of the transfer operator .
in fact , an msm provides 
 an estimate of this transfer operator by estimating a finite number 
 of the eigenfunctions as a sum of indicator functions , where the indicator 
 functions define the state decomposition in phase space . however , 
 the msm is not the only technique that can approximate the transfer 
 operator .
the eigenvalues of the transfer operator correspond 
 to a relaxation time scale , ti , which can be measured via the autocorrelation function of 
 a trajectory { x0,x1, ... ,xt , ... } projected onto 
 the corresponding eigenfunction:7 in recent work from no and 
 nske , it was shown that a variational 
 principle can be derived for the eigenvalues of the transfer operator , 
 meaning that given a mean - centered , integrable function , , the autocorrelation of f at time  is always 
 less than the autocorrelation of the first 
 nonstationary eigenfunction , 1 , of the transfer 
 operator:8 this is powerful because it implies 
 that one can find the slowest 
 dynamical processes in a simulation without building an explicit 
 estimate of the operator itself , but by directly estimating 
 its eigenfunctions .
have successfully constructed estimates of the transfer operator s 
 top eigenfunctions in the span of a prespecified library of basis 
 functions by leveraging this variational bound . here , we take a different 
 approach that does not require a predefined basis set . in
what follows , we are working in the space 
  and will assume that it is a euclidean space with some dimension , d. therefore , if x and y are in 
  , then xy = xy denotes the inner product 
 and x  y = xy the outer product between the two length - d column vectors .
time - structure based independent 
 component analysis ( tica ) is an unsupervised learning method for finding 
 projections that decorrelate slowly . given a markov chain in phase 
 space : { xt   } , 
 tica attempts to find a vector v   such 
 that the autocorrelation of xt projected onto v is maximized ( eq 9).9here ,   and  is some lag time greater than 
 zero .
to compute the expectation values in eq 9 , assume we can sample m pairs of points in phase 
 space , { ( xt  ,yt   ) } , where xt is sampled from the equilibrium 
 distribution , ( ) , and yt from the transition probability density , p(xt, ;  ) . then , 
 we can estimate the time - lag correlation matrix , covariance matrix , 
 and mean vector using a maximum likelihood approach with a multivariate 
 normal approximation ( see si section
let1011where xt = xt , with 
 mean given by12 given these matrices , eq 9 becomes13it can be shown that the solution to eq 13 is 
 also a solution to the generalized eigenvalue 
 problem:14 furthermore , the remaining 
 eigenvectors from eq 14 are the slowest projections 
 constrained to being uncorrelated 
 with the previously found solutions .
this method has proven 
 to be very useful for constructing msms 
 that reproduce the slowest time scales in a simulation .
however , more importantly , since eqs 7 and 9 are the same , the solutions to the tica problem 
 ( tics ) are also solutions to the variational problem and therefore 
 provide estimates of the slowest eigenfunctions of the transfer operator .
these estimates , however , are crude since the 
 tics are constrained to be linear functions of the input coordinates , 
 but if we can extend the tica problem to include nonlinear degrees 
 of freedom then we will be able to estimate the true eigenfunctions 
 more accurately .
the simplest approach 
 to extending the tica problem beyond linear combinations of coordinates 
 in  is to use a mapping function , :  v , which expands the original representation of the data 
 into a higher - dimensional feature space , v. we can 
 then calculate the linear tics in the feature space , v , but since  is arbitrary , these tics can in fact be nonlinear 
 functions of the coordinates in . the dimensionality of v , however , could be enormous ( or perhaps infinite ) , and 
 so , this scheme is impractical except in very simple cases . under 
 certain conditions ,
however , it is possible to perform analysis in v without ever actually calculating the representations , 
 (x ) , by employing what is known as the  kernel 
 trick .  for illustration
, we reproduce a simple example 
 from schlkopf et al . of the kernel 

let   and define the mapping 
 function between 
 phase space , , and feature space , v , 
 to be (x ) = ( x12 , , x22 ) .
because of this definition , 
 it is easy to show that15equation 15 indicates 
 that we can calculate the dot product in v in one 
 of two ways .
we can either compute (x ) and (y ) and explicitly calculate their dot product or we can simply 
 square the dot product of the vectors in .
either route gives us the same result , 
 but the latter has the advantage of being much easier in cases where 
 the feature space s dimension is very large 
imagine 
 using a polynomial of degree much greater than two in eq 15or infinite .
kernels are used commonly 
 in machine learning when extending linear 
 problems to provide nonlinear solutions . any problem that can be rewritten entirely in terms of inner products 

can be extended to use the kernel trick , and in the following , we 
 show how to rewrite the tica problem in this way . remember 
 that the tica problem searches for a projection that maximizes the 
 autocorrelation function of a projection onto a vector v. our goal is to rewrite the optimization problem ( eq 13 ) in terms of only inner products and then apply the kernel 
 trick to solve the tica problem in high - dimensional feature spaces . 

this assumption is only violated when one coordinate can be written 
 as a linear combination of the other coordinates . in these cases , 

we need only project the data into a nonredundant space and continue 
 the analysis there .
we start with a set of m pairs of points in phase space , { ( xt   , yt   ) } , where xt is sampled from the equilibrium distribution , ( ) and yt from the transition probability 
 density , p(xt, ;  ) . consider a function :  v that maps points in phase space ,  , into a higher - dimensional 
 feature space , v. now we want to solve the tica problem 
 in this new space v. the solution , v  v , is a solution to the eigenvalue problem 
 given in eq 14 .
however , the covariance and 
 time - lag correlation matrices are calculated in terms of the vectors 
 in v : here , we have used (x ) 
  (y ) to denote the outer product of the 
 two vectors (x ) , (y )  v. for now , assume that the (xt ) values and (yt ) values are centered ( see supporting information 
 ( si )  centering data in the feature space  )
. it can be shown that the solution , v , is in the span 
 of the xt values and yt values , which means that we 
 can write v as a linear combination of the input data 
 points ( see si 
let   such that16we will rewrite the tica problem in v in terms of a gram matrix k , where k has 2 m rows and 2 m columns , 
 defined as17where 
 (x)(y ) is the inner 
 product between the two vectors (x),(y )  v. in 
 other words , the gram matrix can be broken into four blocks:18such 
 that each block is an m  m matrix and19202122let imm denote the m - dimensional 
 identity matrix and 0mm denote an m by m matrix 
 of all zeros . then , let r be a 2 m  2 m matrix given by eq 23.23 in the si  derivation of the 
 ktica solution ,  we show that the numerator and denominator 
 of the rayleigh quotient in eq 9 can be rewritten 
 in terms of  and k:2425and so the tica problem in the feature space v can 
 be written totally in terms of inner products:26which is itself a rayleigh quotient and so 
 the solutions satisfy the generalized eigenvalue problem:27 however , solving such problems can be challenging because 
 the matrix kk can have zero eigenvalues , which correspond 
 to directions 
 in v that have zero variance associated with them . 

these directions are either due to redundant coordinates ( in v ) or under sampling and so are not desirable solutions 
 to the tica problem . however , these zero variance directions lead 
 to infinite eigenvalues in eq 27 .
this issue 
 is not exclusive to the kernel tica problem and is faced in most kernel 
 learning methods .
several solutions have been reported in the cca , 
 lda , and svm literature that introduce a regularization term .
we use the method applied to kernel - cca in mika et al . but realize that another form of regularization 
 may be even more useful .
the result is a modified objective function 
 that penalizes the solutions according to the l2 norm on .28 this new 
 rayleigh quotient can be written as the solution to a 
 modified generalized eigenvalue problem:29 for large 
 enough  , the result is a positive definite matrix , kk+i2m2 m , which allows us to compute the solutions reliably . 

we refer the reader to other kernel methods that use a similar regularization 
 scheme for a discussion , as well as the results section two - dimensional brownian dynamics . as in the case of linear tica , there are 
 many solutions to the ktica problem ( herein called ktics ) , where each 
 is the slowest projection subject to being uncorrelated with all of 
 the slower ktics .
let v  v denote the i ktic and let   denote the i eigenvector of eq 29 corresponding to v : because of the variational principle , 
 the ktics are direct estimates of the transfer operator s right 
 eigenfunctions : the first eigenfunction ( 0(x ) = 
 1  x   ) is not estimated by ktica 
 as it is stationary , and so the second eigenfunction is estimated 
 by the first ktic .
however , to directly compare to the eigenfunctions 
 of the transfer operator , we need to ensure that the solutions to 
 the ktics are normalized correctly .
the eigenfunctions of the 
 transfer operator are orthonormal with 
 respect to the inner product defined in eq 5 .
since the first eigenfunction is equal to one for all points x , this means that the eigenfunctions have these properties : 1 .
all dynamic eigenfunctions ( i > 0 ) have 
 zero 
 mean ( ) . since the gram matrix is centered in 
 the feature space ( see si section centering 
 data in the feature space ) , the ktics have zero mean .
the 
 eigenfunctions are uncorrelated and are normalized such 
 that they have unit variance ( cov(i(x ) , j(x ) ) = i , j = ij ) . as shown in the si section derivation of the ktica solution , the variance of the i ktic , is given by so we can normalize each  such 
 that the variance is one .
arbitrary points in  can be 
 projected onto the ktics as 
 this operation can be written entirely in terms of inner products:31we emphasize here 
 that it is possible to evaluate 
 the value of v(x)and therefore i+1(x)without ever specifying  or v. in a sense , the ktics behave 
 like a kernel density estimator , where the value at some point x depends on how close that point is to the training points . to solve the ktica 
 problem
, we need to calculate the gram matrix , k , 
 given in eq 18 , which requires the choice of 
 a kernel function that specifies inner products in the feature space . 

although there are a variety of kernel functions available , for this 
 work we chose the commonly used gaussian kernel:32here , x and y correspond 
 to two conformations of a molecule . technically , these vectors are 
 points in phase space and represent the positions and momenta of all 
 atoms in the simulation . however , in most bimolecular simulations , 
 we do not care to describe the rotation and translation of the molecule . 
 moreover , because we only have finite sampling ,
if we attempt to describe 
 rotation and translation in our models , then we will likely do so 
 only by sacrificing accuracy in describing the conformational change 
 that we are interested in . for this reason , for each application we 
 select an internal set of coordinates ( e.g. , all dihedral angles or 
 contact maps of a protein ) to use in eq 32 . 

these coordinates will depend on the particular system being studied , 
 and so we leave their specification to the application sections below . with the kernel given in eq 32 ,
the gram 
 matrix is defined as33 the solutions to eq 29 provide an explicit 
 form of the estimates of each eigenfunction i+1:34 this expression highlights the analogy to a gaussian kernel 
 density 
 estimator , where we can evaluate i+1 at any point x based on the eigenvectors from eq 29 and how close x is to each point in 
 the training set .
in addition to the kernel , we must choose 
 the correlation lag time , 
  , and the regularization strength , . in previous studies , 
 tics were largely unperturbed by the choice of  ; hence , for each of our example systems , we fix 
 a  and ask for the best model given this choice . therefore , 
 there are two free parameters we must specify before we can solve 
 eq 29 :  and . bayes theorem 
 demonstrates that the probability of a model given the data is proportional 
 to the probability of the data given the model , which means that we 
 can select a model based on maximizing the model s likelihood 
 of generating new data . here ,
a  model  is the set of 
 ktics from eq 29 solved with particular values 
 for  and . therefore , our strategy will be to solve 
 eq 29 over a range of values for  and 
  and then pick the model that assigns the highest probability 
 to new data .
as the ktics are estimates of the transfer operator s 
 eigenfunctions , 
 we can define a likelihood based on eq 4 .
let k be a positive integer , then let { ( at  ,bt  )}t  =  1n be pairs of test points 
 that were not used in computing the ktics , where at is sampled from the true equilibrium distribution , 
 ( ) , and bt from 
 the true transition probability density , p(at, ; k ) .
then , the probability of observing all of these transitions can 
 be written in terms of the infinitely many eigenfunctions of the transfer 
 operator : however , we do not have access to infinitely 
 many eigenfunctions , so we instead truncate the expansion at the top m + 1 eigenfunctions to calculate an approximate probability , p.35 this approximation is better when m is large , 
 or when { i : i > m } are very small . since the equilibrium probability is 
 not dependent on either of  and  , we can ignore the 
 ( ) terms when comparing models .
this is especially clear 
 when looking at the approximate log - likelihood:36 for all applications , 
 we calculate a set of ktics for a variety 
 of values of  and  on a subset of the full data set 
 ( the training set ) .
then , eqs 34 and 36 can be used to evaluate the approximate likelihood 
 on { ( at , bt)}t  =  1n  the test set that was initially left out of the calculation . 

the values of  and  that are selected are those that 
 are most likely according to eq 36 .
one 
 technical issue when calculating these probabilities is that 
 since we truncate the infinite expansion at relatively few ktics , 
 the resulting probability is not guaranteed to be positive . in practice
, 
 this can be difficult to handle but can be avoided by only testing 
 pairs of points that are separated by many multiples of . in 
 this way , the faster components we ignore in eq 36 have equilibrated and so the truncated sum is a better approximation 
 of the true ( infinite ) expansion .
admittedly , this problem likely 
 indicates that defining such a probability based solely on the slowest 
 eigenfunctions of the transfer operator is not possible .
there may 
 in fact exist a better scheme for evaluating these models ; however , 
 in the results below , we show that , at least for simple systems , the 
 likelihood approach does work .
computationally , solving the 
 ktica problem requires finding the 
 solution to eq 29 , which is a generalized eigenvalue 
 problem .
however , the matrices on either side are 2 m  2 m , where m is the number of data points in the training set .
this is a computational 
 challenge for all kernel methods , as the solution becomes more difficult 
 to compute as more data is used . here
. however , there are many sophisticated 
 solutions that have been developed for other kernel - learning methods , 
 and these will certainly be useful in the future .
we simulated 
 a markov jump process on a one - dimensional , four - well potential as 
 suggested in prinz et al .
the potential 
 given by37was discretized into 100 points between 1 
 and 1 , and a transition matrix was defined such that at each point , 
 there was nonzero probability to move left , move right , or stay at 
 the current point . in each case , the probability to move from point xi to xj was proportional to exp[(v(xj)v(xi ) ) ] .
from the first 1000 pairs of transitions separated by 10 time steps 
 were used in the ktica method . the euclidean distance between points 

was used in conjunction with the gaussian kernel defined in eq 32 to calculate the gram matrix .
the approximate 
 likelihood scores were calculated for 1000 pairs of points sampled 
 from the second trajectory separated by 10 timesteps . in this simple 
 case ,
there were many models that were roughly the same quality as 
 judged by the likelihood scores ( figure 1 ) .
approximate 
 log likelihood evaluated on test data provided a similar 
 score for most models we tested .
models using a value of 0.01 for 
  or a value of 10 for  were given very low 
 scores , while all other models seemed to be good enough .
the best 
 model according to the likelihood evaluated on test data correctly 
 determined the top three eigenfunctions for this system .
the ktics 
 are shown in blue while the actual eigenfunctions are shown in red . 

it is illustrative to remember that the tica solution is trivial in 
 this one - dimensional case . even in this simple case ,
the nonlinearity 
 in the ktics prove to be essential for estimating the transfer operator s 
 eigenfunctions . since the generating 
 process was itself an msm ,
there are three slow 
 eigenfunctions , corresponding to transitions across the three barriers 
 on this system .
the best model correctly determined the slowest three 
 ktics and are illustrated in figure 1 .
these 
 results illustrate the value of nonlinearity in estimating the eigenfunctions 
 of the transfer operator .
in fact , since the dynamics take place in 
 a one - dimensional space , there is only one solution to the regular 
 tica problem , which only ( crudely ) estimates the slowest eigenfunction 
 of the true transfer operator .
we also investigated the mller 
 potential , and simulated dynamics 
 that were governed bywhere  = 10 , kt = 15 , 
 and r(t ) is a 
 delta - correlated gaussian process with zero mean , and v(x ) was defined aswhere a = ( 1,1,6.5 , 
 0.7 ) ; b = ( 0 , 0 , 11 , 0.6 ) ; c = ( 10,10,6.5 , 
 0.7 ) ; a = ( 200,100,170 , 15 ) ; x = ( 1 , 0,0.5,1 ) ; y = ( 0 , 
 0.5 , 1.5 , 1 ) as suggested by mller and brown . using the euler  maruyama method and a time step 
 of 0.1 , we produced two trajectories of length 10 time 
 steps .
the initial positions were sampled via a uniform distribution 
 over the box : [ 1.5,1.2 ]  [ 0.2,2.0 ] . the 
 ktica algorithm was used to analyze the first trajectory by ( regularly ) 
 sampling 2000 pairs of transitions separated by 10 timesteps .
the 
 cross - validation scores were calculated for 1000 pairs of transitions 
 separated by 50 timesteps using the top two ktics .
as there are essentially 
 three local minima in the mller potential , the top two eigenfunctions 
 should correspond to population transfer between the three wells . 

the model selected by the likelihood criterion correctly estimates 
 these eigenfunctions ( figure 2 ) while numerous 
 other models do not .
ktica models for the mller potential simulations 
 are not 
 as well - behaved as the one - dimensional case .
there were a few models 
 that assigned a nonpositive value for the likelihood of test data 
 ( white regions above ) .
however , the maximum likelihood model correctly 
 identified the top two eigenfunctions , which separate the three major 
 minima in the potential energy surface .
these eigenfunctions are depicted 
 in the bottom two panels , with the red and blue color corresponding 
 to the value of the ktic .
the solutions are highly 
 dependent on the choice of the  , 
 which indicates the necessity of regularization .
if  is too 
 small , then the ktics focus on only a few points that are not sampled 
 often in the data set .
when it is too large , the eigenvalues are underestimated 
 and the eigenfunctions are distorted from the correct functions . the 
 approximate likelihood evaluated on test data chooses an intermediate 
 value of  whose ktics correctly separate the three wells in 
 the mller potential ( figure 3 ) .
three models 
 with  = 0.5 but varying  values are 
 depicted above .
these results illustrate the necessity of using regularization . 
 when  is too small ( top ) ,
the solutions focus on only a few 
 points in the data set that are not sampled heavily , while when it 
 is too large , the solutions are distorted ( bottom ) .
fortunately , the 
 likelihood based model selection criterion picks an intermediate value 
 of  ( middle ) whose ktics are better estimates of the actual 
 eigenfunctions .
a single trajectory of alanine 
 dipeptide ( ala1 ) was simulated using 
 the amber99sb - ildn force field in implicit 
 solvent ( gbsa ) for 500 ns using openmm6.0 with a time step of 2 fs at 300 k using a langevin integrator .
the 
 trajectory was subsampled every 250 fs for analysis and split in half , 
 with the first half serving as the training data and the second serving 
 as the test data .
the ktica calculation was performed by regularly 
 sampling 2000 pairs of points separated by 1.25 ps from the training 
 data .
the gaussian kernel was used in conjunction with a distance 
 metric defined on intra - atomic distances . for each conformation ,
the 
 distance between all pairs of heavy atoms were computed , which formed 
 a vector of length 45 .
then , the distance between two conformations 
 was the euclidean distance between these vectors . for the likelihood 
 calculation
, we regularly sampled 5000 pairs of points separated by 
 6.25 ps to ensure that the approximate likelihoods defined by the 
 top two ktics were non - negative for most models . we should not expect 
 this approximate probability to accurately describe motions at 1.25 
 ps but would hope that as one looks further in time the probability 
 based on the top two eigenfunctions becomes more accurate .
since alanine 
 dipeptide has been the subject of many computational studies , it is 
 known that the slowest two relaxations occur along the axes of the 
 ramachandran plot of the two dihedral angles ,  and . 
 the most likely model recapitulated these two slow processes ( figure 4 ) .
we reanalyzed the two 100 s md trajectories 
 of the fip35 ww 
 domain performed by shaw et al .
the ktica 
 problem was solved by analyzing 10  000 pairs of transitions 
 separated by 50 ns from the first trajectory .
( note that because of 
 the size of the data set , the pairs of transitions overlapped one 
 another .
this is , in a sense , the same as using the sliding window 
 approach for counting transitions when building an msm , which is common 
 practice . )
the gaussian kernel defined in eq 32 was used with a distance metric defined on the protein s 
 contact map .
each conformation was represented as a 528-element vector 
 of intramolecular residue  residue distances . where all pairwise 
 distances were tracked between residues that were at least three amino 
 acids away from each other .
the distance between a pair of residues 
 was computed by taking the distance between the two closest heavy 
 atoms in the respective residues .
the likelihood defined by the slowest 
 ktic alone was evaluated on 1000 pairs of transitions separated by 
 500 ns .
this time was long enough to assign non - negative likelihoods 
 to a majority of the models we tested .
as in the one - dimensional 
 case , there are several models that are 
 very good according to the cross - validation scores , while models with 
 small values of  and  perform poorly .
alanine dipeptide 
 has been studied extensively , and it is known that the slowest two 
 relaxations occur along the axes of the ramachandran plot .
the projection 
 of the test data set into these axes shows that the top two ktics 
 correctly identify the slow motions in the alanine dipeptide simulations .
likelihood defined by only the slowest ktic 
 was evaluated for the 
 ww domain test data .
again , small values of  and  produced 
 likelihoods that were often negative for test data ; however , there 
 is a region in which all likelihoods were positive .
the maximum likelihood 
 model s slowest ktic corresponded to the folding transition , 
 as it was correlated with rmsd to the native state .
the time scale 
 associated with this ktic was 4.4 s , which agrees with previous 
 analysis of this same data set .
the slowest eigenfunctions of 
 the transfer operator for the ww 
 domain are not known , but for the maximum likelihood model , the slowest 
 ktic is physically reasonable as it separated folded and unfolded 
 states in the simulation ( figure 5 ) .
additionally , 
 the ktic had a corresponding time scale of 4.4 s , which is 
 consistent with previous analysis of this data set .
though the tica method has proven to be remarkably successful for 
 building state decompositions for msm analysis , its solutions 
 are limited because they are constrained to be linear functions of 
 the input coordinates .
for this reason , the extension to nonlinear 
 functions by using the kernel trick represents a substantial improvement 
 that allows ktica to provide estimates of the transfer operator s 
 eigenfunctions directly , without the use of an msm .
we have shown here that the method can accurately calculate the 
 slowest eigenprocesses of simple systems but admit that there is still 
 work to be done before it can be routinely used for analyzing arbitrarily 
 complex dynamics .
computationally , the matrices in eq 29 scale with the amount of data used , which means that as the 
 systems become larger and more complicated and more data is needed , 
 the eigenvalue problem will become increasingly difficult to solve . 

however , this problem is not unique to ktica and , in fact , is ingrained 
 in all kernel learning methods . for this reason , we expect significant 
 improvements to be made rapidly as these methods have been extensively 
 studied in machine learning .
 the allure of a molecular dynamics 
 simulation is that , given a 
 sufficiently accurate force field , it can provide an atomic - level 
 view of many interesting phenomena in biology . 
 however , the result 
 of a simulation is a large , high - dimensional time series that is difficult 
 to interpret . 
 recent work has introduced the time - structure based 
 independent components analysis ( tica ) method for analyzing md , which 
 attempts to find the slowest decorrelating linear functions of the 
 molecular coordinates . 
 this method has been used in conjunction with 
 markov state models ( msms ) to provide estimates of the characteristic 
 eigenprocesses contained in a simulation ( e.g. , protein folding , ligand 
 binding ) . here , 
 we extend the tica method using the kernel trick to 
 arrive at nonlinear solutions . 
 this is a substantial improvement as 
 it allows for kernel - tica ( ktica ) to provide estimates of the characteristic 
 eigenprocesses directly without building an msm . 
Introduction
Theory
Results
Conclusions