Analysis of Dynamic Efficacy Profile of Corneal Cross-Linking: Role of Oxygen and Rate Constants

Purpose:To explore (theoretically) the key parameters and their influence on the time profiles of photosensitizer (riboflavin), free radicals, singlet oxygen, oxygen and the efficacy of corneal collagen crosslinking (CXL)in both type-I and oxygen-mediated type-II mechanisms, specially the role of oxygen and the initiator regeneration. Methodology:Coupled kinetic equations are derived and numerically solved under the quasi-steady state condition for the 2-pathway mechanisms of CXL. The key parameters explored include (bI, V, Q', K, K',Q,P) and their influence on the time profiles of photosensitizer (riboflavin, C), radicals (R), singlet oxygen(S), oxygen (X) and efficacy (E), parameters of (K,K',Q) define the relative strength of type-I and type-II process. The oxygen depletion profile, X(t), and the associated singlet oxygen, S(t), depend on the parameters of V, Q' and the initial value of oxygen. The coupling strength given by (bI) governs almost all profiles, where b is an effective absorption parameter and I is the UV light intensity.Results:Our numerical method for CXL dynamic profiles demonstrated the following important features: (i) Type-I and type-II coexit in CXL, in the presence of oxygen. However, there is no type-II when oxygen is depleted or in a condition without oxygen. (ii) Type-I with bimolecular termination, the radical R(t) is a function of [K'(bIgC)]0.5, leading to the steady-state efficacy given by a scaling law of 1/(bI)0.5, in contract to that of type-II which is almost independent to the light intensity. (iii) The depletion rate (2 to 5 minutes) of X(t) is much faster than that of C(t) (10 to 20 minutes), (iv) Thepure type-II profile, has a transition point from straight line to saturating curve and matches the depletion point of singlet oxygen S(t). (v) Improved CXL efficacy of type-I and type-II may be achieved by external supply of photoinitiator (riboflavin) and oxygen, respectively.


INTRODUCTION
Corneal collagen crosslinking (CXL) have been one of the important clinical subjects in ophthalmology. The clinical outcomes of CXL and biomechanical properties of corneal tissue after CXL have been reported and summarized in a book [1]. The safety and efficacy issues of CXL have been reported theoretically [2][3][4][5][6][7][8][9][10], with recent articles reported by Lin et al [11,12]. The critical parameters influencing the efficacy of CXL include: initial concentration and diffusion depth of riboflavin (RF) (for type-I CXL) and oxygen (for type-II CXL), quantum yield, UV light intensity, dose and irradiation duration [11]. Most of the previous models [2][3][4][5][6] are not accurate due to the oversimplified assumptions of constant RF profiles, or non-depleted RF, or UV light intensity following the simple Beer-Lambert law (BLL). Standard (Dresden) protocols were revised for faster (accelerated) CXL based on Bunsen-Roscoe law (BRL) having a limited validation of UV maximum intensity [12]. Controversial efficacy issues of Dresden versus accelerated corneal crosslinking (A-CXL) have been discussed recently by Lin [12] and a concentration-controlled method (CCM) to improve the efficacy of A-CXL was also proposed [10].
Schumacher et al [3] reported the non-oxygen-mediated (NOM) type-I CXL, in contrast to Kling et al [5] claiming that oxygen-mediated (OM) type-II played the critical role of CXL efficacy. Furthermore, Kamaev et al [2] claimed that CXL is NOM-type-I dominant, while the OM-type-II only plays a limited and transient role. If Kling et al [5] were correct, then all the reported results of epi-on CXL and accelerated CXL would not be possible, since only minimum initial oxygen supply is available and the resupply (diffusion) of oxygen takes about 10 minutes [2].
Since the first human data of Wollensak et al in 2003 using the so-called Dresden protocol [1,13], the efficacy of accelerated and standard CXL were reported clinically for the roles of RF concentration and oxygen [14][15][16][17][18][19][20][21][22], where strategy for improved CXL efficacy via enhanced oxygen were reported [17][18][19][20][21][22]. However, due to the lack of reliable protocols, many controversial issues remain as debating, specially the role of oxygen in type-I CXL [12]. Further more, there is no clinical or theoretical systematic studies on the role of the reaction rate constants which govern the efficacy of type-I and type-II CXL.
The present article, for the first time, will explore (theoretically) the key parameters and their influence on the time profiles of photosensitizer (riboflavin), free radicals, singlet oxygen, oxygen and CXL efficacy. The important feature and new finding of CXL based on these numerically produced curves will then be discussed.

MATERIALS AND METHODS
Both type-I and type-II reactions can occur simultaneously, and the ratio between these processes depends on the type of photosensitizers (PS) used, the concentrations of PS, substrate and oxygen, the kinetic rates involved in the process, and the light intensity, dose, PS depletion rate etc. [17]. Detailed kinetic of type-II only, and type-I only was published in my prior work [10] and [11], respectively. Typical depletion time of oxygen is about 5 to 15 seconds, for light intensity of 30 to 3 mW/cm 2 , per measured data of Kamaev et al [2], and takes about 10 minutes for the oxygen to be resupplied or replenished to about 1/3 of its initial state.
As shown in Figs. 1, the CXL process is described as follows [17]. The ground state RF molecules (C) are excited by the UV light to its excited state (C*) and triplet excited state (T). In type-I process, (T) could interact directly with the stroma collagen substrate (M) for crosslinking, and produces a radical (R) and regenerate the initiator (C). T could also interact with the ground state oxygen, [O2], to form reactive superoxide anion radicals [O-] (not shown in Fig.1). For type-II process, T interacts with [O2] to form oxygen singlet [ 1 O2], which could be relaxed to its ground state oxygen [O2], or crosslink the stroma substrate (M). Besides, CXL, oxygen singlet could be used to kill bacteria for the treatment of corneal karatatitis or for anti-cancers.
The kinetic equations (based on the kinetic chart of Fig. 1) for the concentration of various components are shown were shown previously by Lin et al [ ]. We will use the following short-hand notations: C and T for the RF ground and excited triplet state; R for the active radical, S for the singlet oxygen    where b=a'q83.6  , with q=k2/(k1+k2) is the quantum yield of T; a' is the extinction coefficients of RF ;  being the UV light wavelength. k2 is the relaxation rate of radical (S). Re is a regeneration term. V and Q' are the coupling constants for R and X, T and X, respectively. K, K' and Q are the crosslink constant for T, R and S. The steady-state values are given by: T=bIgC, S= g'Q' TX; with g=1/(V+ KM+ k4X); g'=1/(k6+ k1C+ QM); and the radical (R) is more complex given by the solution of [ ] where G= VX+ K'M and H= KMT; with T=bIgC. Solving for R, we obtain Analytic formulas of R is available under two special cases.
Case (i) for unimolecular termination dominant, or G>>k'H, we obtain R= K(bIgCM/G)(1-0.5H/G), which is a linear increasing function of bIgC/G, for first-order with 0.5H<<G. In this case, there is a oxygen inhibition (OIH) effect which reduces the radical (R) and the efficacy, because G is an increasing function of oxygen (or X), G= k"X+ K'M. Case (ii) for bimolecular termination dominant, with H>> GR, we obtain, R=[H/k'] 0.5 . a nonlinear function of [K(bIgC)M] 0.5 , a square root function. In contrast to case (i), the OIH plays no role in case (ii) We note that g=1/(k"+ KM+ k4X) defining the lifetime of T (given by 1/g); Eq. (2) includes an oxygen source term given by P=(1-X/X0) P0, with a maximum rate constant P0, where (1-X/X0) is included to avoid the negative value of oxygen. We note, in Eq. (1), -bIC is the RF depletion, which is compensated by a regeneration term, Re, such that dC/dt=-bIC(1-Re)= -(k1CS -k'RX)=0, in the absence of oxygen, or X=S=0. This was the conventionally believed situation that there is no RF depletion in type-II pathway. In fact, in a pure type-I case, with X=S=0, the perfect compensation (with Re=1) is always valid, but not for type-II case. For more complex schemes, this perfect cycle might not be met [18]. High efficacy requires a long lifetime of R and T (or small g). The conversion eq. (3) includes both terms for type-I (KT and K'R) and type-II (QS).
We note that Eq. (1) is much more complex and complete than that of Kamave et al [2], which is our special case when C(t) is a constant (using a continuing resupply of RF), k1C=0 in g', and V=0 in Eq. (2). Our Eq. (3) has 3 terms for crosslink, but Kamave et al [2] ignored the K'R, and assumed monomer conversion is only due to the coupling of T and M, and the coupling of singlet oxygen (S) and M. Kinetic Equations of Schumacher et al [3] and Semchishen et al [4] are limited to the type-I conversion, K'R, and ignored the oxygen-mediated term, QS in our Eq. (3). They also ignored the Re cycle effects. The modeling of Kling [5] is based on Kamave [2], but only showed the algorithm for numerical calculations, there is no analytic formulas. In addition, most previous models, except that of Lin et al [ ], assumed a constant light intensity under a constant C(t). In general, light intensity is an increasing function of time due to the depletion of C(t). Comparing to the above described previous modeling, our modeling, shown by Eq. (1) to (3), is the most complete and accurate one.
The dynamic light intensity is given by [7,11] where, a=83.6  ,  being the UV light wavelength; a'=204 (1/%/cm) and b' are the extinction coefficients of RF and the photolysis product, respectively; Q=13.9 (1/cm) is the absorption coefficient of the stroma at the UV wavelength. We note that the initial light intensity is also given by its deep dependence (z) by I0(z)=I0 exp(-bI0C'z), with C'(t) is an averaged C(z,t) over z. In the following discussion, we will focus on the dynamic profiles on the corneal surface, or for z=0. The z-profiles are much more complex and require numerical simulations per published by Lin et al [ ].

RESULTS AND DISCUSSION
In the following, we will show numerical results based on Eq. (1) to (3) with steady-state radical R given by Eq. (4). The key parameters to be explored include (bI, Q', K, K',Q,P,V) and their influence on the time profiles of photosensitizer (riboflavin, C), radicals (R), singlet oxygen(S), oxygen (X) and efficacy (E). The important feature and new finding of CXL based on these numerically produced curves will then be discussed. We note that, as shown in Eq. (3), the parameters of (K,K',Q) define the relative strength of type-I and type-II process. The oxygen depletion profile, X(t), and the associated singlet oxygen, S(t), depend on the parameters of V, Q' and its initial value of oxygen. The coupling strength given by (bI) governs almost all profiles, where b is an effective absorption parameter and I is the UV light intensity. The parameters of T, R and S are all increasing function of (bI) which serves as a collective effect, i.e. for the same product value of (bI) leads to same efficacy and dynamic profiles. For examples the pair value of (b, I) = (0.1,5), (0.2, 2.5), (0.5, 1.0), and (5.0, 0.1) having the same product value of 0.5 and therefore lead to the same impact, The role of each of the 5 parameters (bI, Q', K, K', Q) will be explored by varying only one, having the others fixed.
We note that the solution of Eq. (3) leads to the efficacy (E) defined by E=1 -M(t)/M(0), with M(0) being the initial monomer (or collagen matrix) concentration. The overall efficacy has three components given by the time integral of KT, K'R (for type-I) and QS (for type-II). Therefore, setting of K'>0 with Q=0 defines pure type-I efficacy, whereas K'=0 and Q>0 defines pure type-II; and mixed process defined by K', Q' >0. shows that both initiator, C(t), and oxygen, X(t) are decreasing function of (bI), as indicated by Eq. (1) and (2). However, the oxygen X(t) depletion (within 2 to 4 minutes) is much faster than that of initiator C(t), within 15 to 20 minutes. These features are in consistent with measured data [ ] that higher light intensity leads to faster depletion of both C and X, for a given absorption constant b. Fig. 3. shows he time profile of free radical R(t) and singlet oxygen S(t),in which higher (bI) leads to lager R(t) and S(t) and larger area covered by the curve, or time integral, which defines the efficacy (E), as shown by Eq. (3), Fig. 4 also shows the saturation feature of the efficacy profile resulted form the complete depletion of C(t), R(t) and S(t) at steady-state. We note that for the case of pure type-I (with Q=X=0), the solution of Eq. (3) is given by the time integral of the factor [bIC0exp(-bIt)] 0.5 , where the 0.5 power leads to its time integration to a factor of 2[1-C0exp(-0.5bIt)]/(bI) 0.5 , which has a steady-state proportional to 1/(bI) 0.5 , or a decreasing function of light intensity. These uncommon feature is also shown by our numerical curves shown in Fig. 4, that higher light intensity has a faster transient value but a lower steady-state efficacy. This feature is consistent with the clinical data that accelerated CXL (with intensity >18 mW/cm 2 ) is less efficient than standard low power CXL (at 3 mW/cm 2 ). We note that this feature only occurs in type-I case (given by R), whereas pure type-II (give by QS, with K'R=0) has the similar steady-state efficacy and almost independent to the light intensity, as shown by Fig. 5, in which there is no crossing curve like Fig. 4. Mathematically, our previous analytic formulas also showed that radical R (for type-I) is proportional to the square-root of (bI), whereas singlet oxygen S (for type-II) is linearly proportional to (bI), therefore their time integral (or the efficacy) have very different feature, as shown by Fig. 4 and 5.     Fig. 7 shows the the role of Q for the case of pure type-I case in the absence of oxygen (with Q=0), and in the presence of oxygen with various conversion rate of S, Q=0.5 and 4.0, in Eq. (3). The efficacy profiles E(t) without oxygen (with Q=0, the red-curve) shows a smooth exponential profiles; as also shown by Fig.6; whereas the profiles for combined type-I and type-II (with Q>0), show the mixed of straight line (for the transient stage) and curved after the depletion of oxygen (at about 3 minutes). Mathematically, this transition feature could be realized by the solution of Eq. (3) that E(t) is slowdown after oxygen (X) and singlet oxygen (S) are depleted, (or when X=S=0), as also shown by profiles of R(t) and S(t) in Fig.3. We note that the efficacy is insensitive to the values of V and Q', therefore we did not show the profiles for various V and Q'.  [20,21]. Fig. 9 shows the efficacy time profile for various initial oxygen X0 with fixed bI=(left figure); and various bI=(0.005, 0.01, 0.02) with fixed X0=1.0. Fig. 10 shows the time profile of oxygen X(t) and efficacy E(t) for various external oxygen supply with P=(0, 0.05, 0.1).

The role of regeneration
In the previous cases, we have used the non-regeneration case with Re=0, leading to a maximum depletion without the cyclic effects for the initiator. A perfect regeneration case (when Re=1.0, dC/dt=0, presenting a constant initiator and leading to highest efficacy. The degree of regeneration can be explored by Re=(0,0.5, 0.8) with the associate profiles of C(t), radical R(t) and efficacy E(t) are shown in Fig. 11, in which higher regeneration leads to higher C,R and E. Fig. 11 also shows the efficacy E(t) for pure type-I case (with Q=0), and mixed case (with Q=1.0) which has higher efficacy than Q=0.

Summary of up-dated CXL Features
From the numerical produced curves shown by Fig. 2 to Fig. 11, the key parameters (bI, Q', K, K',Q,P) and their influence on the time profiles of photosensitizer (riboflavin, C), radicals (R), singlet oxygen(S), oxygen (X) and efficacy (E) are explored. As shown in Eq. (3), the parameters of (K,K',Q) define the relative strength of type-I and type-II process. The oxygen depletion profile, X(t), and the associated singlet oxygen, S(t), depend on the parameters of V, Q' and its initial value of oxygen. The coupling strength given by (bI) governs almost all profiles, where b is an effective absorption parameter and I is the UV light intensity. The parameters of T, R and S are all increasing function of (bI) which serves as a collective effect, i.e. for the same product value of (bI) leads to same efficacy and dynamic profiles.
From our numerical analysis, the important features of type-I and type-II CXL and each of the key parameters of (bI, Q', K, K',Q,P) are summarized and compared as follows: (a) Type-I and type-II coexit in CXL, in the presence of oxygen. However, there is no type-II when oxygen is depleted or in a condition without oxygen. (b) Type-I with bimolecular termination, the radical R(t) is a nonlinear square-root function of [K'(bIgC)] 0.5 , which leads to a unique feature that the steady-state efficacy is a decreasing function of light intensity and follows a scaling law of 1/(bI) 0.5 , in contract to that of type-II which is almost independent to the light intensity; comparing Fig. 4 and 5. (c) The depletion rate (2 to 5 minutes) of oxygen X(t) is much faster than that of the initiator C(t), (10 to 20 minutes), as shown by Fig. 2; in which the associated radicals R(t) and S(t) have profiles cutoff at the time when C(t) and X(t) are completely depleted, as shown by Fig. 3. (d) Efficacy profile E(t) is given by the time integral of the crosslink radicals. Fig. 4 shows the profiles for pure type-I case (with Q=0), in which the areas covered by R(t) defines the steady-state of E(t) and inverse proportional to (bI) 0.5 , but its rising rate of is proportional to (bI). (e) Fig. 5 shows the pure type-II profile (with K=0), the profile has a transition point (TP) from straight line to saturating curve given by the depletion point of singlet oxygen S(t); where TP=(2,2.5,4.2) minutes matching the potion of S(at t=TP)=0. (f) Profiles of Fig. 6 demonstrate the role of the rate constants K and K' given in Eq. (3). In contrast to Fig. 4 having steady-state inverse proportional to (bI), it is proportional to K and K'. Similarly, Fig. 7 shows the role of rate constant Q. Furthermore the profiles transition points (at about 2 minutes is defined by the depletion point of R(t) and S(t). (g) Fig. 8 shows the role oxygen X(t) having steady-state value nonlinearly proportional to its initial value X0 (shown by right figure). We note that oxygen is required for oxygen-mediated type-II but it is not required in in type-I. Therefore, type-II only plays a limited and transient state role for t<t0, with t0 being the depletion time of oxygen. The conventionally believed concept that oxygen is a must-element for CXL which is controversial as discussed by Lin [12]. (h) Fig. 9 compares the roles of X0 and (bI) on the oxygen profiles, which also show the depletion points associate with the TP of Fig. 8. (i) It was known that depletion of photoinitiator C(t) could be compensated by C(t) self-regeneration, or by external supply during the CXL [10]. For the case of perfect regeneration case, or Re=1 in Eq. (1), C(t) is a constant due to the catalytic cycle. Fig. 11 shows that slower depletion of C(t) leads to larger area covered by R(t) which defines the efficacy E(t). Fig. 11 also compares the efficacy for pure type-I and mixed process with enhanced efficacy due to type-II. (j) It was clinically reported that both type-I and type-II coexist In the transient stage until the oxygen is depleted; then type-I dominates before the oxygen is resupplied or replenished [2]. Improved efficacy was also demonstrated clinically by external supply of oxygen [ ]. Fig. 10 demonstrates this feature theoretically for pure type-II case (with K'=0). As predicted and shown by Fig. 2, the depletion of C(t) is much slower than that of oxygen. Therefore, at the time oxygen is depleted, (or type-II reaches its steady-state efficacy), approximately 60% to 80% of initiator is still available for type-I process. We further note that (not shown here) in a pure type-II CXL, in the absence of oxygen supply (or P=0), higher light intensity has a faster rising curve, but all intensities reach the same steady state value. However, for P>0, higher intensity has lower steady state value due to the faster oxygen depletion-profiles, shown in our previous model [19,20].
We finally note that Enhanced efficacy could be achieved by various strategies including external supply of oxygen and/or initiator, increased light intensity (but limited to transient state), and higher rate constants via for example, co initiators or additives [23]. Furthermore, many of our theoretical new findings require further clinical support using our predicted features as the protocol guidance.

CONCLUSION
Our new numerical simulation showed that oxygen (and singlet oxygen) play a limited and transient role in the process, in consistent with that of Kamave [2]. In contrary, Kling et al [3] believed that type-II is the predominant mechanism, which however, conflicting with the epi-on CXL results. For both type-I and type-II, a transient state conversion (crosslink) efficacy in an increasing function of light intensity (or dose), whereas, its steady-state efficacy (in pure type-I case) is a deceasing function of light intensity governed by a scaling law of 1/(bI) 0.5 . The depletion rate (2 to 5 minutes) of oxygen X(t) is much faster than that of the initiator C(t) (10 to 20 minutes). The pure type-II profile, has a transition point (TP) from straight line to saturating curve and matches the depletion point of singlet oxygen S(t). Enhanced efficacy could be achieved by various strategies including external supply of oxygen and/or initiator, increased light intensity (but limited to transient state), and higher rate constants via for example, co initiators or additives.

CONSENT
It is not applicable.

ETHICAL APPROVAL
It is not applicable.