Open AccessArticle

Charge Recombination Kinetics of Bacterial Photosynthetic Reaction Centres Reconstituted in Liposomes: Deterministic Versus Stochastic Approach

Emiliano Altamura

Paola Albanese

Pasquale Stano

Massimo Trotta

Francesco Milano

^4,*

and

Fabio Mavelli

^1,*

Department of Chemistry, Università degli studi di Bari “Aldo Moro”, Via Orabona 4, 70125 Bari, Italy

Department of Biological and Environmental Sciences and Technologies (DiSTeBA), University of Salento, Strada Provinciale Lecce-Monteroni, Ecotekne, 73100 Lecce, Italy

IPCF-CNR Istituto per i Processi Chimico Fisici, via Orabona 4, 70125 Bari, Italy

⁴

Institute of Science of Food Production (ISPA), National Research Council of Italy (CNR), Strada Provinciale Lecce-Monteroni, Ecotekne, 73100 Lecce, Italy

Authors to whom correspondence should be addressed.

Data 2020, 5(2), 53; https://doi.org/10.3390/data5020053

Submission received: 17 May 2020 / Revised: 8 June 2020 / Accepted: 9 June 2020 / Published: 12 June 2020

Download

Browse Figures

Versions Notes

Abstract

In this theoretical work, we analyse the kinetics of charge recombination reaction after a light excitation of the Reaction Centres extracted from the photosynthetic bacterium Rhodobacter sphaeroides and reconstituted in small unilamellar phospholipid vesicles. Due to the compartmentalized nature of liposomes, vesicles may exhibit a random distribution of both ubiquinone molecules and the Reaction Centre protein complexes that can produce significant differences on the local concentrations from the average expected values. Moreover, since the amount of reacting species is very low in compartmentalized lipid systems the stochastic approach is more suitable to unveil deviations of the average time behaviour of vesicles from the deterministic time evolution.

Keywords:

reaction centres; liposomes; charge recombination; ODE; stochastic simulations

1. Introduction

The bacterial photosynthetic Reaction Centre (RC) is a transmembrane pigment-protein complex [1,2,3], that, upon photoexcitation, generates a charge-separated state among a bacteriochlorophyll dimer: the primary donor (D) and two ubiquinone-10 complexes (Q) located in two protein sites indicated as Q_A (the primary acceptor) and Q_B (the secondary acceptor), respectively [4]. The electron is reversibly exchanged between the two ubiquinones with a thermodynamic constant, named L_AB, whose magnitude is related to the stability of Q⁻ in the site Q_B respect to the stability in the site Q_A. In presence of the exogenous electron donor, a reduced cytochrome c₂ molecule, a second electron reaches the ubiquinone in the Q_B site, leading to its full reduction to ubiquinol, upon the uptake of two protons from the external aqueous solution. The ubiquinol then leaves the protein from the Q_B site and it is replaced by a new ubiquinone molecule of the membrane pool [5]. Both ubiquinone and ubiquinol reversibly bind and unbind to the Q_B site of the protein [6], with binding constants, named K_Q and K_QH2, respectively. It is worthwhile to remark, that the ubiquinone in the Q_A site is more strongly bound to RC and it cannot leave this site without inactivating the protein in the physiological conformation. In isolated RC, in absence of any exogenous electron donor, if the charge-separated state is generated by a short saturating light irradiation, it is possible to follow in the dark the charge recombination (CR) reaction at 865 nm, where the reduced D has an absorption peak that is bleached upon photo-oxidation to D⁺. The CR kinetics can be used to retrieve kinetic and thermodynamic information about the energetics of the electron acceptor complex and the ubiquinone binding. Being both RC and ubiquinone strongly hydrophobic, the system is typically studied in detergent micelles, where one RC protein interacts with zero, one or few ubiquinone-10 molecules [7]. This situation, however, is far from the physiological condition, where 11 dimeric and 2 monomeric LH1-RC complexes [8] are present in the intracytoplasmic membrane (ICM) together with a different amount of ubiquinones. An environment that better mimics the physiological membrane is represented by liposomes whose lipid composition can be tailored to closely resemble those of ICMs, although not all the proteins and molecules constituting the photosynthetic apparatus of the bacterium are present [9,10,11]. It is important to remark that, for this reason, in the less crowded membrane of liposomes, the ubiquinone molecular diffusion can be different from what is observed in the bacterial ICM. In a previous work [10], RCs were reconstituted in lecithin vesicles (mainly constituted by the zwitterionic palmitoyl-oleoyl-phosphatidylcholine-POPC) at high Q/RC ratio. L_AB and K_Q were determined accounting for the ubiquinone concentration polydispersity among different liposomes by building a ubiquinone concentration distribution function: P(N_Q), that was used to fit the experimental CR curves. In a subsequent work [11], RCs were reconstituted in vesicles made of either POPC or the negatively charged palmitoyl-oleoyl-phosphatidylglycerol (POPG) at variable Q/RC ratio. In this case, L_AB and K_Q were determined using a deterministic approach by building a system of ordinary differential equations (ODE), describing all reactions that occur in the system, assuming a uniform ubiquinone concentration.

In the present theoretical work, we try to overcome the limits of the deterministic kinetic description of compartmentalized reacting systems [12] by using a stochastic approach [13,14,15,16] to take into account and evaluate both intrinsic and extrinsic random effects. Intrinsic random effects are mainly due to the fact that reacting molecules are not all in the same energetic conditions and this can induce fluctuations in the occurring time of reactions. These fluctuations influence the time course of the reaction systems the lower the number of molecules that react [16]. Therefore, systems confined in space, such as liposomes that encapsulate enzymes [17], could be highly affected by intrinsic stochastic effects. On the other hand, the polydispersity of a vesicle suspension, both in sizes and compositions, is a source of the so-called extrinsic random effects [18]. Examples reported in recent literature show how the preparation of lipid vesicles can produce a distribution of entrapped molecules that can affect significantly the time behaviour of chemically reacting liposomes [19,20].

Following the deterministic approach and using the experimental POPC decay curves already published by our group [11], this work focuses the attention on how stochastic effects can influence the average time behaviour of a vesicle population determining significant displacements from the deterministic CR time course. This will be done performing Monte Carlo simulations that allow to calculate the CR reaction curves averaging the stochastic behaviour of a large number of single vesicles, using a software platform [21] that implements the Gillespie’s stochastic simulation algorithm [22,23]. Different examples of Monte Carlo methods for the simulation of stochastic kinetics of chemically reacting systems have been also reported in earlier works [24,25].

2. Methods

2.1. The Kinetic Mechanism

Figure 1 reports the assumed kinetic mechanism describing the charge recombination reaction occurring when isolated RC is hit by a saturating excitation light radiation [6].

The mechanism describes the spontaneous charge recombination processes that take place in the dark, i.e., when the light is switched off: t = 0. All the steps reported in the mechanism are considered elementary, i.e., single event steps. For this reason, the deterministic reaction rates can be expressed by means of the mass action law, given by the product of the kinetic constant times all reactant concentrations raised to a power of their own stoichiometric coefficient. In the bottom line, the reversible ubiquinone association/dissociation to/from the

Q_{B}

site in the dark (neutral state) is reported. On the other hand, the ubiquinone is assumed always present and blocked in the Q_A site when the protein complex is photoactive in its physiological conformation. When the system is hit by a suitable saturating light radiation, the species in the dark equilibrium are instantaneously brought in the charge-separated state and all the possible electron transfer, charge recombination and association/dissociation steps take place when the light is switched off. Indeed, the charge separation generated by the light irradiation occurs within a picosecond range and it can be considered instantaneous respect to all the backward reactions of the RC mechanism. Indeed, the photochemical excitation is not explicitly studied in this work.

2.2. Experimental Outcomes and Kinetic Rate Constant Estimation

The time evolution of the CR reactions can be monitored following the different absorbance (

Δ A

) between neutral and charge-separated species. Therefore, the normalized absorbance change

Δ A / Δ A_{0}

is proportional to the molar fraction of all the species in the charge separated state χ(D⁺):

χ (D^{+}) = χ (D^{+} Q_{A}^{-}) + χ (D^{+} Q_{A}^{-} Q_{B}) + χ (D^{+} Q_{A} Q_{B}^{-})

The experimental curves are usually fitted by a bi-exponential phenomenological decay function:

χ(D⁺) = A_F exp(−k_F t) + A_S exp(−k_S t) + c

(1)

that is the sum of a fast decay process (A_F exp(−k_F t)), and a slower charge recombination path (A_S exp(−k_S t)), plus a constant term c accounting for background noise (typically < 2% of the overall signal).

Indeed, if the Q_B pocket is empty or inhibited by herbicides, such as terbutryn, CR occurs only by direct tunnelling from

Q_{A}^{-}

with a mono-exponential decay. In this case, the time constant k_F is quite independent from the RC surroundings and approximates the direct charge recombination constant (k_F = k_AD ≈ 10 s⁻¹), A_F ≈ 1.0 and the slow decay contribution is quite negligible A_S ≈ 0.0.

On the other hand, at saturating-Q_B site ubiquinone concentrations, the decay is still mono-exponential, but slower, with the slow time constant approximately equal to: k_S = k_AD/(1 + L_AB) (where L_AB = k_AB/k_BA), and A_S ≈ 1.0 while A_F ≈ 0.0. In this condition, the k_S value is influenced by the protein solubilizing environment [26,27], since the membrane lipid composition affects the association/dissociation processes of Q at the Q_B site. In fact, the mentioned retardation arises because the direct route from

D^{+} Q_{A} Q_{B}^{-}

is so slow that can be negligible [26] and the population of

D^{+} Q_{A}^{-} Q_{B}

is lowered by the fast reversible electron transfer to

D^{+} Q_{A} Q_{B}^{-}

When the ubiquinone concentration is in the sub-saturating-Q_B site range, the CR reactions become generally biphasic. The first exponential accounts for the direct charge recombination (fast decay) and its time constant (k_F) is independent from the occupancy of the

Q_{B}

-site. The time constant k_S of the second exponential (slow decay) results from a complex kinetic constant combination of electron transfers and ubiquinone binding (

k_{in}^{*}

k_{in}

) and unbinding (

k_{out}^{*}

k_{out}

) from the

Q_{B}

-site [27]. Since in these conditions the different solubilizing environment can mainly affect the ubiquinone release kinetic constant k_out, it is possible to define different ubiquinone exchange regimes: the slow one (k_out ≪ k_AD), typical of the detergent environment [28], and the fast one (k_out ≫ k_AD), typical of the reverse micelle environment [29]. When k_out ≈ k_AD the kinetics of CR reactions shows a simultaneous

Q_{B}

-site population dependence of both A_S and k_S, like in the case investigated in this work [11,30]. It is important to remark that, in these experimental conditions, only the optimization of theoretical curves obtained by numerical integrations of the complete mechanism can extract the kinetic constant values from the mixed-order kinetics [27].

The previously described best-fit procedure has been applied to the study of CR reactions of Q_B-depleted RCs reconstituted in zwitterionic (POPC) liposomes adding increasing amounts of ubiquinone [11]. The charge recombination experimental traces were recorded at 865 nm and the normalized absorbance decay has been fitted by using the bi-exponential function Equation (1). According to a new method introduced [11], K_Q = k_in/k_out (calculated considering all concentrations referred to the bilayer phase) was found to be 1.42 × 10³ M⁻¹ for RCs embedded in POPC liposomes at 25 °C. This allows to estimate the upper values of k_in by the equation k_in = f_c·ξ·4 π r·d_Q·N_A where f_c is the RC-Q collision frequency, ξ is the Q_B site entrance fraction of RC surface, r = 30 Å is the radius of the protein, N_A is the Avogadro’s number, and d_Q is the ubiquinone diffusion coefficient in the vesicle bilayer [11]. This results in k_in ≤ 2.18 × 10⁵ M⁻¹ s⁻¹ and k_out ≤ 154 s⁻¹ obtained from k_in/K_Q. The value of k_AD is fixed at those obtained for Q_B-depleted RCs (k_F = 10 s⁻¹) and the k_BD value is fixed at 0.06 s⁻¹ [26]. Finally, the value of k_AB for POPC at 25 °C is experimentally available [31], resulting k_AB = 8.15 × 10³ s⁻¹, while k_BA is calculated to be 574 s⁻¹ from the relation k_AB/L_AB with L_AB = 14.2 obtained from the value of k_S. All these values are reported in Table 1 as guess values.

2.3. The Kinetic Ordinary Differential Equation Set

The set of ordinary differential equations (ODE) associated to the mechanism of charge recombination reported in Figure 1 is the following:

{\begin{matrix} {dx / dt = - (k}_{AD} {+ k}_{in} {q) x + k *}_{out} y \\ {dy / dt = k}_{in} {qx + k}_{BA} {z - (k}_{AD} {+ k}_{AB} {+ k *}_{out}) y \\ {dw / dt = k}_{out} {u + k}_{AD} {x - k}_{in} qw \\ {dz / dt = k}_{AB} {y - (k}_{BA} {+ k}_{BD}) z \\ {du / dt = k}_{in} {qw + k}_{AD} {y + k}_{BD} {z - k}_{out} u \\ {dq / dt = - k}_{in} {q (w + x) + k *}_{out} {y + k}_{out} u \end{matrix}

(2)

Here

, {k'}_{AD}

and

k_{AD}

have been assumed equal because the presence of a ubiquinone in the Q_B site does not perturb the charge recombination rate from Q_A⁻. Indeed, the distance of the Q_B site from the Q_A site is very large. The assumption

k_{in}^{*}

k_{in}

, instead, is only based on practical reasons. Indeed, the pairs

k_{in}^{*}

k_{out}^{*}

and

k_{in}

k_{out}

are highly correlated, since their ratios give the values of K*_Q and K_Q, respectively. Moreover,

k_{in}^{*}

and

k_{in}

exhibit the highest values being second order kinetic constants of fast processes and this makes them less sensitive and very hard to optimize simultaneously, bringing always the best-fit procedure towards a complete set of kinetic parameters with a physically poor meaning. Consequently, we decided to optimize separately only

k_{out}

and

k_{out}^{*}

and use these parameters to determine the ubiquinone exchange regime in the neutral and in the illuminated state, respectively.

The species concentrations have been made dimensionless with the following transformations:

\begin{array}{l} x = [D^{+} Q_{A}^{-}] / [RC]; \\ y = [D^{+} Q_{A}^{-} Q_{B}] / [RC]; \\ z = [D^{+} Q_{A} Q_{B}^{-}] / [RC]; \\ w = [{DQ}_{A}] / [RC]; \\ u = [{DQ}_{A} Q_{B}] / [RC]; \\ q = [Q]_{free} / [Q]_{T}; \end{array}

The concentrations reported in square brackets are calculated as local concentrations, i.e., molar concentrations referred to the volume of the liposomal lipid membrane, for all the species involved in the charge recombination mechanism. Indeed, the RC transmembrane proteins and ubiquinone molecules are present only in the bilayer of liposomes, since they are highly hydrophobic species. The local concentration [X] can be obtained dividing the average number of X molecules per liposome

\bar{N_{X}}

by the Avogadro’s number N_A and by the volume of the lipid membrane

V_{M e m} = \frac{4}{3} π [R^{3} - {(R - τ)}^{3}]

, R being the liposome hydrodynamic radius and τ the bilayer thickness (4 nm). The average number

\bar{N_{X}}

can be estimated as the ratio between the bulk concentrations

c_{X}^{T}

divided by the concentration of liposomes

c_{V e s}

that can be determined by using the bulk concentration of phospholipid

c_{L}^{T}

and the phospholipid polar head area

α_{L}

. Therefore, local concentrations can be calculated with the following formula:

[X] = \frac{\bar{N_{X}}}{N_{A} V_{M e m}} = \frac{c_{X}^{T}}{N_{A} V_{M e m} c_{V e s}} = \frac{c_{X}^{T}}{c_{L}^{T}} \frac{{8 π R}^{2}}{α_{L} N_{A} V_{M e m}} = \frac{c_{X}^{T}}{c_{L}^{T}} \frac{{6 R}^{2}}{N_{A} α_{L} [R^{3} - {(R - τ)}^{3}]}

(3)

where 8 π R² is the total area of the interphase between the bilayer membrane and the aqueous solution. By using the hydrodynamic radius of liposomes of 55 ± 6 nm obtained by the DLS analysis, and considering bulk concentrations

c_{R C}^{T}

= 1.0 µM of RC and

c_{L}^{T}

= 1.0 mM of POPC, respectively, the local RC concentration results [RC] = 1.2 ± 0.3 mM. For each liposome, the local concentration [Q]_T of the total ubiquinone was calculated assuming that the bulk ratio Q/RC is kept the same at least on average in each compartment at the end of the preparation procedure:

c_{Q}^{T}

c_{R C}^{T}

= [Q]_T/[RC] [32].

It is important to remind that the overall ubiquinone concentration [Q]_T is given by the sum of the free ubiquinone [Q] and of the RC with the Q_B site occupied: [Q]_T = [Q] +

[{DQ}_{A} Q_{B}]

[D^{+} Q_{A}^{-} Q_{B}]

[D^{+} Q_{A} Q_{B}^{-}]

, because the ubiquinone in the Q_A site cannot leave the protein pocket.

2.4. Numerical Integrations

The integration of the ODE set in Equation (2) has been performed by a homemade code written in Matlab ver.2019b by using a numerical procedure suitable for stiff problems: ode23s, based on a modified Rosenbrock formula of order 2 [33].

The initial concentration values of the reacting species are determined by the pre-flash equilibrium concentrations of DQ_A, DQ_AQ_B and Q according to the dark binding equilibrium constant K_Q = k_in/k_out and setting the [RC] = 1.0 × 10⁻³ M and [Q]_T = [RC] ×

c_{Q}^{T} / c_{R C}^{T}

. Moreover, light irradiation produces the charge separated state, i.e., [D⁺Q_A⁻]₀ = [DQ_A]_Eq, and [D⁺Q_A⁻Q_B]₀ = [DQ_AQ_B]_Eq on a timescale much faster than CR decay. All the other species are null. The local initial concentration values for the species different from zero are reported in Table 2.

Finally, the CR decay against time t is obtained from the ODE set solutions as the sum of all species concentrations in the charge separation state divided by [RC], the local concentration of the all reaction centres present in the vesicle membrane:

χ(D⁺)_t = ([D⁺Q_A⁻]_t + [D⁺Q_A⁻Q_B]_t + [D⁺Q_AQ_B⁻]_t)/[RC].

2.5. Optimization Procedure

For data-fitting-problems, the Matlab procedure lsqcurvefit has been used by setting the Levenberg-Marquardt algorithm [34,35]. Ad hoc Matlab codes have been developed for the optimization of the numerical curves and the experimental data by using the Matlab’s toolboxes. In the Supplementary Information, the code for the numerical integration is reported as an example.

2.6. Stochastic Simulations

The time evolution of a population of RC-embedding liposomes after the light radiation switched off has been simulated by using the ENVIROMENT program: a software platform developed in the past years to study the stochastic kinetics of compartmentalized, self-reproducing, chemically reacting systems [21]. This platform, based on the Stochastic Simulation Algorithm introduced by Gillespie in the seventies [22,23] and improved in following years [36], is suitable to run Monte Carlo time simulations of a population of reacting compartments dispersed in size and composition [37,38,39]. According to the Gillespie algorithm, the time evolutions obtained by simulations can be considered as sub-sample of the statistical ensemble, rigorously described by the Kinetic Master Equation [13,14,15]. Therefore, they can be analysed to get insights on the average time behaviour of the liposome population along with the standard deviations, estimated as the square root of the average quadratic displacements from the mean.

It is worthwhile to remark that the stochastic approach is based on the Mean Field and on the Continuous well-Stirred Tank Reactor (CSTR) assumptions. The first one considers all the reacting molecules being in the same physical conditions at least in average. Because of this, the average stochastic time behaviour of the reacting systems must converge to the deterministic one in the thermodynamic limit, i.e., when the number of the reacting molecules tends to infinite [40]. In order to satisfy this requirement, the stochastic probability coefficient p_ρ must be related to the deterministic kinetic rate constant k_ρ of the elementary reaction ρ by the following equation [13,14,40]:

p_{ρ} = \frac{k_{ρ}}{{(V_{ρ} N_{A})}^{μ_{ρ} - 1}}

where µ_ρ is the molecularity of the ρ-th elementary reaction. The CSTR states that the reacting system is chemically homogeneous, since it is continuously stirred. For real reacting systems, this means that the diffusion processes are much faster than chemical reactions, a requirement that is usually satisfied in normal operative conditions, i.e., at room temperature and in aqueous solutions [41].

The stochastic simulations performed with populations of vesicles of different sizes have been carried out always assuming perfectly spherical and monodisperse vesicle populations, i.e., with the same hydrodynamic radius R. Therefore, the contribution to the displacements from the population average time trends of the reacting species concentrations due to the dimensional dispersion of the vesicles was not considered. Instead, in this work, a comparison among the average time behaviours of different dimensional classes of spherical vesicles with different radius has been performed.

On the other hand, to evaluate the role of the extrinsic stochastic effects, the vesicle populations have been filled at the beginning of each simulation run assuming a uniform distribution U(N_x) or a random distribution, i.e., approximately Gaussian N(n|N_x,

\sqrt{N_{x}}

)dn and Poisson P(n|N_x) probabilities. In the case of the uniform distribution, the local concentration [X] of each species has been converted into an integer number of molecules N_X encapsulated in the membrane volume V_memb by the formula N_X = [X] N_AV_memb, being N_A the Avogadro’s number. This molecular number was set in all of the vesicles of the population for the species X. In the case of a random distribution, fixed M as the number of vesicles, the total number of X molecules to be entrapped in the vesicle population has been calculated: N^T_X = M[X]N_AV_memb, and then these molecules have been distributed by chance among all the M vesicles drawing pseudo-random numbers uniformly distributed in the range 1-M.

In Figure 2 a comparison between simulated distribution of the initial number of molecules in the vesicle population and the expected theoretical probabilities is reported in the case: [RC]= × 10⁻³ M and

c_{Q}^{T} / c_{R C}^{T} = 0.1

In the first case, all the vesicles of the considered population start from the same initial conditions and the observed differences between the deterministic and simulated trends are due only to the intrinsic stochastic fluctuations of the different reactive processes. In the second case, the initial distribution of the solutes is carried out considering deviations from the average value due to completely random fluctuations and, therefore, with a standard deviation equal to the square root of the average number of molecules. This approach allows to consider the uncertainty of the initial conditions of the different vesicles, i.e., the extrinsic effects, that can differentiate the deterministic behaviour from the stochastic one.

3. Results and Discussion

According to a method introduced in a previous paper [11], the kinetic constants of the complete mechanism were estimated from the bi-exponential best fit values obtained for different ubiquinone concentrations. These kinetic data are reported in Table 1 as guess values, used in the present work, as input data for the next optimization step. In this step, the numerical solutions of the complete ODE set are used for determining all the rate constants of the mechanism in Figure 1 as optimization values.

In Figure 3, the comparison between the experimental data (black lines) and the optimized theoretical curves (red lines) is reported, while the kinetic constant best fit values are listed in Table 1.

Although the agreement between ODE solutions and experimental data is quite satisfactory, the residual analysis shows that statistically significant displacements between experimental data and numerical curves are exhibited by all the four studied cases.

To account for these displacements, stochastic simulations have been then carried out. The aim of this study is to understand in which conditions a significant difference between stochastic and deterministic behaviour could emerge, due to intrinsic and extrinsic random effects. This can rationalize the difference observed between the experimental data and the deterministic curves.

Simulations have been, then, performed for different POPC/RC/Q ratios and different sizes of perfectly spherical and dimensionally monodisperse vesicle populations, by assuming both uniform and random distribution of RC and ubiquinone molecules. The average concentration trends over time were then calculated and compared with the deterministic numerical solutions of the ODE system reported in Equation (2). Both stochastic simulations and deterministic curves have been determined by using the kinetic constant best-fit values reported in Table 1.

3.1. Stochastic Simulations of a Single Vesicle with Increasing Membrane Volume

In all the following simulations, the local concentration of Q_B-depleted RC has been set to [RC] = 1 × 10⁻³ M, while the reaction volume and the ratio

c_{Q}^{T}

c_{R C}^{T}

have been changed, as shown below:

\begin{matrix} \begin{array}{l} Analysed \\ systems \end{array} & {\begin{array}{l} V_{40} = 7 {. 265 \times 10}^{- 20} d m^{3}, \times 1, \times 10, \times 10^{2}, \times 10^{3} \\ \frac{c_{Q}^{T}}{c_{R C}^{T}} = 0.1, 0 . 3, 0 . 5, 0 . 7, 1 . 0, 3 . 0 \\ [R C] = 0.001 M \end{array} \end{matrix}

In particular, the reactive volume starts from V₄₀, i.e., the volume of 4.0 nm thick membrane of a lipid vesicle with a hydrodynamic radius R = 40 nm, and increases being multiplied by the factor f ranging from 1 to 10³. Therefore, the obtained values correspond to the bilayer volumes of hypothetical spherical vesicles of radius 40, 120, 380 and 1200 nm respectively. In the lowest V₄₀, N_RC is 44 and it increases according to the vesicle membrane volume.

In Figure 4, the comparisons between the simulated trends and the solutions of the ODE system for the decay of the photo excited species χ* against time are reported for different membrane volumes and different

c_{Q}^{T}

c_{R C}^{T}

ratios.

The figure shows that, for all the

c_{Q}^{T}

c_{R C}^{T}

ratios studied, the simulated trends converge towards the deterministic curves as the reaction volume increases. This is confirmed by the calculated values of the square root of the mean square displacement (RMSD) between the ODE and stochastic curves (STO):

R M S D = \sqrt{\frac{\sum_{i} {(χ_{i, S T O}^{*} - χ_{i, O D E}^{*})}^{2}}{N}}

(4)

In the V₄₀ case, the stochastic trend is very strongly influenced by random fluctuations that become almost negligible at higher volumes. These results confirm that the stochastic system behaviour tends to the deterministic one in the thermodynamic limits, as expected. Indeed, the RMSD value drops below 1% already for vesicles of R = 122 nm, as shown by the reported values in the legend of Figure 4 (see curves STO × 10¹). However, since 40 nm is around the experimental radius determined by DLS analysis, these graphs clearly show that the time behaviour of real vesicles can be greatly influenced by the intrinsic stochastic effects.

3.2. Stochastic Simulations of Vesicle Populations with Constant Radius and Uniform Solute Distribution

Having shown that intrinsic random fluctuations play an important role, especially for small reaction vesicles where high displacements from the deterministic behaviour can be observed, we investigate if the stochastic time trends tend to the deterministic curves. The simulations are carried out over a 500-vesicle population for different

c_{Q}^{T}

c_{R C}^{T}

values, assuming a uniform solute molecules distribution, to avoid the contribution of extrinsic stochastic effects. All the scanned conditions are listed, as follows:

\begin{matrix} \begin{array}{l} Analysed \\ systems \end{array} & {\begin{array}{l} V_{40} = 7 {. 265 \times 10}^{- 20} d m^{3} \\ \frac{c_{Q}^{T}}{c_{R C}^{T}} = 0.1, 0 . 3, 0 . 5, 0 . 7, 1 . 0, 3 . 0 \\ [RC] = 0.001 M \end{array} \end{matrix}

Figure 5 shows the comparisons between the stochastic average trends (black lines) and the deterministic curves (red lines) calculated by solving the ODE set for all the ratios

c_{Q}^{T}

c_{R C}^{T}

. A 500-vesicle population seems to be large enough to reduce RMSD, in all cases, to values lower than those observed experimentally. Like in the case of the single vesicle simulations (Figure 4), higher displacements are exhibited by higher

c_{Q}^{T}

c_{R C}^{T}

values.

It is worthwhile to remark that the deterministic curves, reported in Figure 5, have been calculated by using initial concentrations determined as simulated local concentrations: [X] = N_X/(N_AV₄₀), i.e., taking into account the rounding effect due to integer numbers of reacting molecules per vesicle. In fact, deterministic curves calculated using strictly the bulk initial concentrations show a worst agreement with the stochastic averages (data not shown). This highlights how, for small reacting compartments, the ODE system solution is very sensitive to uncertainties on the initial concentrations that, therefore, should be determined experimentally, with high accuracy, whenever possible.

3.3. Stochastic Simulations of Vesicles Populations with Gaussian Solute Distribution

The last set of simulations is carried out on a population of 500 monodisperse, spherical vesicles for different radius and different

c_{Q}^{T}

c_{R C}^{T}

ratios, assuming a random Gaussian solute distribution to verify how differences on the local composition can influence the average stochastic behaviour:

\begin{matrix} \begin{array}{l} Analysed \\ systems \end{array} & {\begin{array}{l} R = 20, 30, 40, 50 and 200 nm \\ \frac{c_{Q}^{T}}{c_{R C}^{T}} = 0.1, 0 . 3, 0 . 5, 0 . 7, 1 . 0, 3 . 0 \\ [R C] = 0.001 M \end{array} \end{matrix}

Figure 6A shows the case R = 200 nm and ratio

c_{Q}^{T}

c_{R C}^{T}

= 1.0 in which there is a good agreement between stochastic averages and deterministic curves, with RMSD less than 1% and the confidence bands are very close to the red trace. In Figure 6B,C the cases R = 40 nm and R = 20 nm and the same ratio

c_{Q}^{T}

c_{R C}^{T}

= 1.0 are respectively reported. In both these cases, the discrepancy between the stochastic averages and the deterministic trends is quite evident, despite the fact that the population of vesicles taken into consideration is statistically significant. The lower radius case (R = 20) has the highest RMSD of 4%, a value that is found for all

c_{Q}^{T}

c_{R C}^{T}

ratios (data not shown). Furthermore, the confidence bands became wider as R decreases, pointing out the important effect of fluctuations even for large vesicle populations.

The marked difference observed between the deterministic and stochastic behaviours is essentially due to the possibility of finding liposomes with [Q]_T/[RC] local ratios significantly different from the average value

c_{Q}^{T}

c_{R C}^{T}

. Moreover, these deviations are stronger in vesicle populations with smaller radius. In fact, these liposomes exhibit completely different CR decays (for instance, mono- against bi-exponential decays) that the average on the population cannot compose.

The simulations carried out have shown that in the thermodynamic limit the two approaches return convergent trends for all the species involved in the charge recombination process of the photo-excited reaction centres. This is the case of vesicles with R > 200 nm and of small vesicles (R ≤ 40 nm) if a uniform distribution of reagent species is assumed. However, if a random distribution of solutes is introduced, spherical vesicles with radius in the range of 10–100 nm can exhibit a stochastic behaviour that differs significantly from the deterministic one. This can be explained by the presence of vesicles with local ratios [Q]_T/[RC] significantly different from the average expected value

c_{Q}^{T}

c_{R C}^{T}

4. Conclusions

In this article, the charge recombination kinetics of photo-excited reaction centres reconstituted in lipid vesicles have been analysed with two complementary approaches. The first one, the deterministic approach, is based on the numerical solution of the ODE set associated to the CR kinetic mechanism. It has been applied to the case of RCs embedded into lipid membranes of phosphatidylcholine liposomes to determine, as best-fit values, the kinetic constants of each reactive step in the CR mechanism. The very good agreement, obtained between the optimized ODE solution curves and the experimental data, has been due to the good initial kinetic parameters used by the optimization procedure. Indeed, the starting kinetic parameters have been derived by analysing the best-fit parameters of a bi-exponential phenomenological function and by values reported in literature [11].

Based on the set of the optimized kinetic constants, the CR time evolution of RCs embedded into membranes of a liposome population has been simulated by a Monte Carlo program [21] based on the Gillespie’s Stochastic Simulation Algorithm [22,23,36]. The simulation outcomes have allowed to elucidate the influence of intrinsic and extrinsic random effects on the average time behaviour of the liposome population, due to both the low numbers of RC proteins and Q cofactors embedded in the lipid membranes and the polydispersity of vesicles in terms of composition and size. It is important to remark that, since no assumption on the vesicle size distribution has been done, in this work only the influence of the different liposome size on the random fluctuations have been shown and discussed. As expected, we found that the simulated time evolutions tend, in average, to the deterministic solutions when a population of larger size vesicles is considered and uniform concentrations of embedded RCs and Qs are used. On the other hand, significant differences have been found if random solute distributions are applied for vesicles with R < 40 nm. This can be ascribed to the fact that, in compartmentalized systems with low number of reacting molecules, i.e., in smaller vesicles, the local ratio

c_{Q}^{T}

c_{R C}^{T}

can be, by chance, highly different from the expected macroscopic value. Larger displacement from the average trend can be observed for

0 < c_{Q}^{T}

c_{R C}^{T} < 3

and for smaller vesicles. This can drive towards a completely different time behaviour from what is deterministically predicted in highly size polydispersevesicle suspensions. This effect needs more investigation to be better elucidate and it will be the topic of future studies. Moreover, this kind of studies can be better grounded if experimental data on size and composition distributions of the reacting lipid compartments are available, but unfortunately these information are very hard to get experimentally.

In conclusion, stochastic simulations can be suitable to get insights in the time behaviour of small compartmentalized chemically reacting systems with a randomly disperse composition and highly useful in the characterization of artificial lipid systems that mimic the behaviour of real cells [42,43].

Supplementary Materials

The following are available online at https://www.mdpi.com/2306-5729/5/2/53/s1, supplementary information: Matlab Code for the numerical integration of the ODE set.

Author Contributions

Conceptualization, M.T., F.M. (Francesco Milano) and F.M. (Fabio Mavelli); Data curation, E.A. and F.M. (Fabio Mavelli); Investigation, E.A., P.A. and F.M. (Francesco Milano); Methodology, F.M. (Fabio Mavelli); Software, F.M. (Fabio Mavelli); Supervision, F.M. (Fabio Mavelli); Writing—original draft, F.M. (Francesco Milano) and F.M. (Fabio Mavelli); Writing—review & editing, E.A. and P.S. All authors have read and agreed to the published version of the manuscript.

Conflicts of Interest

The authors declare no conflicts of interest.

References

Bahatyrova, S.; Frese, R.N.; Siebert, C.A.; Olsen, J.D.; van der Werf, K.O.; van Grondelle, R.; Niederman, R.A.; Bullough, P.A.; Otto, C.; Hunter, C.N. The native architecture of a photosynthetic membrane. Nature 2004, 430, 1058–1062. [Google Scholar] [CrossRef]
Allen, J.P.; Williams, J.C. Photosynthetic reaction centers. FEBS Lett. 1998, 438, 5–9. [Google Scholar] [CrossRef] [Green Version]
Allen, J.P.; Feher, G.; Yeates, T.O.; Komiya, H.; Rees, D.C. Structure of the reaction center from Rhodobacter sphaeroides R-26: Protein-cofactor (quinones and Fe2+) interactions. Proc. Natl. Acad. Sci. USA 1988, 85, 8487–8491. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Stowell, M.H.; McPhillips, T.M.; Rees, D.C.; Soltis, S.M.; Abresch, E.; Feher, G. Light-induced structural changes in photosynthetic reaction center: Implications for mechanism of electron-proton transfer. Science 1997, 276, 812–816. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Mcpherson, P.H.; Okamura, M.Y.; Feher, G. Electron transfer from the reaction center of Rhodobacter sphaeroides to the quinone pool - doubly reduced Q_B leaves the reaction center. Biochim. Biophys. Acta 1990, 1016, 289–292. [Google Scholar] [CrossRef]
Crofts, A.R.; Meinhardt, S.W.; Jones, K.R.; Snozzi, M. The role of the quinone pool in the cyclic electron-transfer chain of Rhodopseudomonas sphaeroides: A modified q-cycle mechanism. Biochim. Biophys. Acta 1983, 723, 202–218. [Google Scholar] [CrossRef] [Green Version]
Altamura, E.; Fiorentino, R.; Palazzo, G.; Stano, P.; Mavelli, F. First moves towards photoautotrophic synthetic cells: In vitro study of photosynthetic reaction centre and cytochrome bc1 complex interactions. Biophys. Chem. 2017, 229, 46–56. [Google Scholar] [CrossRef] [PubMed]
Singhaory, A.; Maffeo, C.; Delgado-Magnero, K.H.; Swainsbury, D.J.K.; Sener, M.; Kleinekathöfer, U.; Vant, J.W.; Nguyen, J.; Hitchcock, A.; Isralewitz, B.; et al. Atoms to Phenotypes: Molecular Design Principles of Cellular Energy Metabolism. Cell 2019, 179, 1098–1111. [Google Scholar] [CrossRef] [PubMed]
Altamura, E.; Milano, F.; Tangorra, R.R.; Trotta, M.; Omar, O.H.; Stano, P.; Mavelli, F. Highly oriented photosynthetic reaction centers generate a proton gradient in synthetic protocells the article. Proc. Natl. Acad. Sci. USA 2017, 114, 3837–3842. [Google Scholar] [CrossRef] [Green Version]
Palazzo, G.; Mallardi, A.; Giustini, M.; Berti, D.; Venturoli, G. Cumulant analysis of charge recombination kinetics in bacterial reaction centers reconstituted into lipid vesicles. Biophys. J. 2000, 79, 1171–1179. [Google Scholar] [CrossRef] [Green Version]
Mavelli, F.; Trotta, M.; Ciriaco, F.; Agostiano, A.; Giotta, L.; Italiano, F.; Milano, F. The binding of quinone to the photosynthetic reaction centers: Kinetics and thermodynamics of reactions occurring at the Q_B-site in zwitterionic and anionic liposomes. Eur. Biophys. J. 2014, 43, 301–315. [Google Scholar] [CrossRef] [PubMed]
Mavelli, F.; Altamura, E.; Cassidei, L.; Stano, P. Recent theoretical approaches to minimal artificial cells. Entropy 2014, 16, 2488–2511. [Google Scholar] [CrossRef]
Van Kampen, N.G. Stochastic Processes in Physics and Chemistry, 1st ed.; Elsevier: Amsterdam, The Netherlands, 1981; pp. 166–186. [Google Scholar]
Gillespie, D.T. Markov Processes: An Introduction for Physical Scientists; Academic Press: San Diego, CA, USA, 1992; pp. 332–335. [Google Scholar]
Érdi, P.; Lente, G. Stochastic Chemical Kinetics—Theory and (Mostly) Systems Biological Applications; Springer: New York, NY, USA, 2014. [Google Scholar]
Mavelli, F.; Piotto, S. Stochastic simulations of homogeneous chemically reacting systems. J. Mol. Struct. Theochem 2006, 771, 55–64. [Google Scholar] [CrossRef]
Küchler, A.; Yoshimoto, M.; Luginbühl, S.; Mavelli, F.; Walde, P. Enzymatic Reactions in Confined Environments. Nat. Nanotechnol. 2016, 11, 409–420. [Google Scholar] [CrossRef]
Altamura, E.; Carrara, P.; D’Angelo, F.; Mavelli, F.; Stano, P. Extrinsic stochastic factors (solute partition) in gene expression inside lipid vesicles and lipid-stabilized water-in-oil droplets: A review. Synth. Biol. 2018, 3. [Google Scholar] [CrossRef] [Green Version]
D’Aguanno, E.; Altamura, E.; Mavelli, F.; Fahr, A.; Stano, P.; Luisi, P.L. Physical routes to primitive cells: An experimental model based on the spontaneous entrapment of enzymes inside Micrometer-Sized liposomes. Life 2015, 5, 969–996. [Google Scholar] [CrossRef] [Green Version]
Mavelli, F.; Stano, P. Experiments and numerical modelling on the capture and concentration of transcription-translation machinery inside vesicles. Artif. Life 2015, 21. [Google Scholar] [CrossRef]
Mavelli, F.; Ruiz-Mirazo, K. ENVIRONMENT: A computational platform to stochastically simulate reacting and self-reproducing compartments. Phys. Biol. 2010, 3, 36002. [Google Scholar] [CrossRef]
Gillespie, D.T. A general method for numerically simulating the stochastic time evolution of coupled chemical reactions. J. Comput. Phys. 1976, 22, 403–434. [Google Scholar] [CrossRef]
Gillespie, D.T. Exact stochastic simulation of coupled chemical reactions. J. Phys. Chem. 1977, 81, 2340–2361. [Google Scholar] [CrossRef]
Sipos, T.; Tóth, J.; Érdi, P. Stochastic simulation of complex chemical reactions by digital computer. I. The model. React. Kinet. Catal. Lett. 1974, 1, 113–117. [Google Scholar] [CrossRef]
Sipos, T.; Tóth, J.; Érdi, P. Stochastic simulation of complex chemical reactions by digital computer. II. Applications. React. Kinet. Catal. Lett. 1974, 1, 209–213. [Google Scholar] [CrossRef]
Kleinfeld, D.; Okamura, M.Y.; Feher, G. Electron transfer in reaction centers of Rhodopseudomonas sphaeroides. I. Determination of the charge recombination pathway of D⁺Q_AQ⁻_B and free energy and kinetic relations between Q−AQB and QAQ−B. Biochim. Biophys. Acta 1984, 766, 126–140. [Google Scholar] [CrossRef]
Shinkarev, V.P.; Wraight, C.A. Electron and proton transfer in the acceptor quinone complex of reaction centres of phototrophic bacteria. In The Photosynthetic Reaction Center; Deisenhofer, J., Norris, J.R., Eds.; Academic Press: Cambridge, MA, USA, 1993; pp. 193–255. [Google Scholar]
Agostiano, A.; Milano, F.; Trotta, M. Investigation on the detergent role in the function of secondary quinone in bacterial reaction centers. Eur. J. Biochem. 1999, 262, 358–364. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Mallardi, A.; Palazzo, G.; Venturoli, G. Binding of ubiquinone to photosynthetic reaction centers: Determination of enthalpy and entropy changes in reverse micelles. J. Phys. Chem. B 1997, 101, 7850–7857. [Google Scholar] [CrossRef]
Milano, F.; Agostiano, A.; Mavelli, F.; Trotta, M. Kinetics of the quinone binding reaction at the Q(B) site of reaction centers from the purple bacteria Rhodobacter sphaeroides reconstituted in liposomes. Eur. J. Biochem. 2003, 270, 4595–4605. [Google Scholar] [CrossRef]
Nagy, L.; Milano, F.; Dorogi, M.; Agostiano, A.; Laczko, G.; Szebenyi, K.; Varo, G.; Trotta, M.; Maroti, P. Protein/lipid interaction in the bacterial photosynthetic reaction center: Phosphatidylcholine and phosphatidylglycerol modify the free energy levels of the quinones. Biochemistry 2004, 43, 12913–12923. [Google Scholar] [CrossRef]
Milano, F.; Italiano, F.; Agostiano, A.; Trotta, M. Characterisation of RC-proteoliposomes at different RC/lipid ratios. Photosynth. Res. 2009, 100, 107–112. [Google Scholar] [CrossRef]
Shampine, L.F.; Reichelt, M.W. The MATLAB ODE Suite. SIAM J. Sci. Comput. 1997, 18. [Google Scholar] [CrossRef] [Green Version]
Levenberg, K. A Method for the Solution of Certain Problems in Least Squares. Q. Appl. Math. 1944, 2, 164–168. [Google Scholar] [CrossRef] [Green Version]
Marquardt, D. An Algorithm for Least-Squares Estimation of Nonlinear Parameters. SIAM J. Appl. Math. 1963, 11, 431–441. [Google Scholar] [CrossRef]
Li, H.; Cao, Y.; Petzold, L.R.; Gillespie, D.T. Algorithms and software for stochastic simulation of biochemical reacting systems. Biotechnol. Prog. 2008, 24, 56–61. [Google Scholar] [CrossRef] [PubMed]
Mavelli, F.; Ruiz-Mirazo, K. Stochastic simulations of minimal self-reproducing cellular systems. Phil. Trans. R. Soc. Lond. Ser. B Biol. Sci. 2007, 362, 1789–1802. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Mavelli, F. Stochastic simulations of minimal cells: The Ribocell model. BMC Bioinform. 2012, 13. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Mavelli, F.; Stano, P. Kinetic models for autopoietic chemical systems: Role of fluctuations in homeostatic regime. Phys. Biol. 2010, 7. [Google Scholar] [CrossRef] [PubMed]
Kurtz, T.G. The Relationship between Stochastic and Deterministic Models for Chemical Reactions. J. Chem. Phys. 1972, 57, 2976–2978. [Google Scholar] [CrossRef]
Lente, G.; Ditrói, T. Stochastic Kinetic Analysis of the Frank Model. Stochastic Approach to Flow-Through Reactors. J. Phys. Chem. B 2009, 113, 7237–7242. [Google Scholar] [CrossRef]
Stano, P.; Altamura, E.; Mavelli, F. Novel directions in molecular systems design: The case of light-transducing synthetic cells. Commun. Integr. Biol. 2017, 10, e1365993. [Google Scholar] [CrossRef]
Altamura, E.; Albanese, P.; Marotta, R.; Milano, F.; Fiore, M.; Trotta, M.; Stano, P.; Mavelli, F. Light-driven ATP production promotes mRNA biosynthesis inside hybrid multi-compartment artificial protocells. bioRxiv 2020. [Google Scholar] [CrossRef]

Figure 1. Set of elementary reactions for the dark relaxation of photosynthetic Reaction Centres embedded in liposome membranes.

k_{AD}

and

{k'}_{AD}

are the kinetic constants of the electron transfer from

Q_{A}^{-}

to D⁺ when the

Q_{B}

pocket is empty or occupied, respectively. k_BD is the kinetic constant of the direct charge recombination from

Q_{B}^{-}

to D⁺.

k_{in}^{*} {/ k}_{out}^{*}

and

k_{in} / k_{out}

represent the kinetic constants of the ubiquinone association/dissociation to/from the

Q_{B}

site in the charge separated and neutral state, respectively; k_AB is the kinetic constant for the electron transfer reaction from

Q_{A}^{-}

Q_{B}

, while k_BA is the kinetic constant for the backward reaction.

Figure 1. Set of elementary reactions for the dark relaxation of photosynthetic Reaction Centres embedded in liposome membranes.

k_{AD}

and

{k'}_{AD}

are the kinetic constants of the electron transfer from

Q_{A}^{-}

to D⁺ when the

Q_{B}

pocket is empty or occupied, respectively. k_BD is the kinetic constant of the direct charge recombination from

Q_{B}^{-}

to D⁺.

k_{in}^{*} {/ k}_{out}^{*}

and

k_{in} / k_{out}

represent the kinetic constants of the ubiquinone association/dissociation to/from the

Q_{B}

site in the charge separated and neutral state, respectively; k_AB is the kinetic constant for the electron transfer reaction from

Q_{A}^{-}

Q_{B}

, while k_BA is the kinetic constant for the backward reaction.

Figure 2. Random distribution of the reacting molecules D⁺Q_A⁻, D⁺Q_A⁻Q_B and Q at the beginning of a simulation run by setting [RC] = 1.0 × 10⁻³ M and

c_{Q}^{T} / c_{R C}^{T} = 0.1

. Comparisons between the simulate distributions (bars) and the theoretical probabilities, Gaussian N(n|N_x,

\sqrt{N_{x}}

)dn (blue curves) for D⁺Q_A⁻ and Poisson P(n|N_x) (red points and lines) for D⁺Q_A⁻Q_B and Q, are reported.

Figure 2. Random distribution of the reacting molecules D⁺Q_A⁻, D⁺Q_A⁻Q_B and Q at the beginning of a simulation run by setting [RC] = 1.0 × 10⁻³ M and

c_{Q}^{T} / c_{R C}^{T} = 0.1

. Comparisons between the simulate distributions (bars) and the theoretical probabilities, Gaussian N(n|N_x,

\sqrt{N_{x}}

)dn (blue curves) for D⁺Q_A⁻ and Poisson P(n|N_x) (red points and lines) for D⁺Q_A⁻Q_B and Q, are reported.

Figure 3. Experimental charge recombination traces at 865 nm for POPC at 25 °C: normalized absorbance χ(D⁺) =

Δ A / Δ A_{0}

vs time (black lines) along with the optimized solutions (red lines) of the ODE set reported in Equation (2). Q/RC ratios are: 0, 0.5, 0.8 and 20 going from the fastest decay to the slowest. Initial concentrations have been set, as described in Numerical Integration Section. The best-fit kinetic constants are reported in Table 1. The optimized kinetic parameters reduce the average RMSD of the four curves of about 50%.

Figure 3. Experimental charge recombination traces at 865 nm for POPC at 25 °C: normalized absorbance χ(D⁺) =

Δ A / Δ A_{0}

Figure 4. Time evolution of the photo excited species: comparisons between the single runs of stochastic simulations (black and gray curves) and the numerical solution of the ODE set (red dashed curves). Plots are referred to the same RC concentration: c^T_RC = 1.0 × 10⁻³ M, but to six different values of the ratio c^T_Q/c^T_RC. In each plot simulation outcomes were performed for a vesicle with the same membrane composition, but with an increasing membrane volume, starting from the membrane volume of a 40 nm radius liposome V₄₀ =

7 {. 3 \times 10}^{- 20} d m^{3},

(black curve) and enlarging the membrane volume by 10 to 10³ times. Legends show the RMSD values for the individual stochastic traces calculated with Equation (4).

7 {. 3 \times 10}^{- 20} d m^{3},

(black curve) and enlarging the membrane volume by 10 to 10³ times. Legends show the RMSD values for the individual stochastic traces calculated with Equation (4).

Figure 5. Time evolution of the photo excited species: comparisons between the averages of stochastic simulations (black line with gray confidence bands) and the deterministic curves (red lines) for different

c_{Q}^{T}

c_{R C}^{T}

. The averages were calculated over the stochastic simulation outcomes of a 500 monodispersed vesicle population having a 40 nm radius and a membrane volume V₄₀ =

7 {. 3 \times 10}^{- 20} d m^{3},

the same RC concentration: c^T_RC = 1.0 × 10⁻³ M, but six different values of the ratio c^T_Q/c^T_RC. The RC proteins and Q molecules were distributed uniformly over the vesicle populations.

c_{Q}^{T}

c_{R C}^{T}

. The averages were calculated over the stochastic simulation outcomes of a 500 monodispersed vesicle population having a 40 nm radius and a membrane volume V₄₀ =

7 {. 3 \times 10}^{- 20} d m^{3},

the same RC concentration: c^T_RC = 1.0 × 10⁻³ M, but six different values of the ratio c^T_Q/c^T_RC. The RC proteins and Q molecules were distributed uniformly over the vesicle populations.

Figure 6. Time evolution of the photo excited species: comparison between the averages of stochastic simulations (black lines with grey error bands) and the deterministic curves (red lines) for a population of 500 monodisperse vesicles of 200 nm (A), 40 nm (B), and 20 nm (C) radius, respectively, with the same membrane composition:

c_{R C}^{T}

= 1.0 × 10⁻³ M and

c_{Q}^{T}

c_{R C}^{T}

= 1. The RC proteins and Q molecules were distributed randomly over the vesicle populations according to a Gaussian density/Poisson probability.

c_{R C}^{T}

= 1.0 × 10⁻³ M and

c_{Q}^{T}

c_{R C}^{T}

= 1. The RC proteins and Q molecules were distributed randomly over the vesicle populations according to a Gaussian density/Poisson probability.

Table 1. Kinetic rate constants obtained by the analysis of phenomenological bi-exponential optimized curves (guess values) and by the optimization of the numerical solution curves (best-fit values).

Kinetic Rate Constants
Parameters	Guess Values	Best-Fit Values
k*_in = k_in	≤ 2.18 × 10⁵ M⁻¹ s⁻¹	2.45 × 10⁵ M⁻¹ s⁻¹
k*_out	≤ 154 s⁻¹	114 s⁻¹
k_out	≤ 154 s⁻¹	170 s⁻¹
k_AD = k’_AD	10 s⁻¹	9.7 s⁻¹
k_BD	6.0 × 10⁻²	2.18 × 10⁻²
k_AB	8.15 × 10³ s⁻¹	1.0 × 10⁴
k_BA	574 s⁻¹	544 s⁻¹

Table 2. Initial local concentration for the CRreagents different from zero immediately after the light irradiation.

$c_{Q}^{T} / c_{R C}^{T}$	[D⁺Q_A⁻]₀ = [DQ_A]_Eq	[D⁺Q_A⁻Q_B]₀ = [DQ_AQ_B]_Eq	[Q]
0.1	9.319 × 10⁻⁴	6.809 × 10⁻⁵	3.191 × 10⁻⁵
0.3	8.055 × 10⁻⁴	1.945 × 10⁻⁴	1.055 × 10⁻⁴
0.5	6.932 × 10⁻⁴	3.068 × 10⁻⁴	1.932 × 10⁻⁴
0.7	5.960 × 10⁻⁴	4.040 × 10⁻⁴	2.960 × 10⁻⁴
1.0	4.776 × 10⁻⁴	5.224 × 10⁻⁴	4.776 × 10⁻⁴
3.0	1.677 × 10⁻⁴	8.323 × 10⁻⁴	2.168 × 10⁻³

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

Altamura, E.; Albanese, P.; Stano, P.; Trotta, M.; Milano, F.; Mavelli, F. Charge Recombination Kinetics of Bacterial Photosynthetic Reaction Centres Reconstituted in Liposomes: Deterministic Versus Stochastic Approach. Data 2020, 5, 53. https://doi.org/10.3390/data5020053

AMA Style

Altamura E, Albanese P, Stano P, Trotta M, Milano F, Mavelli F. Charge Recombination Kinetics of Bacterial Photosynthetic Reaction Centres Reconstituted in Liposomes: Deterministic Versus Stochastic Approach. Data. 2020; 5(2):53. https://doi.org/10.3390/data5020053

Chicago/Turabian Style

Altamura, Emiliano, Paola Albanese, Pasquale Stano, Massimo Trotta, Francesco Milano, and Fabio Mavelli. 2020. "Charge Recombination Kinetics of Bacterial Photosynthetic Reaction Centres Reconstituted in Liposomes: Deterministic Versus Stochastic Approach" Data 5, no. 2: 53. https://doi.org/10.3390/data5020053

APA Style

Altamura, E., Albanese, P., Stano, P., Trotta, M., Milano, F., & Mavelli, F. (2020). Charge Recombination Kinetics of Bacterial Photosynthetic Reaction Centres Reconstituted in Liposomes: Deterministic Versus Stochastic Approach. Data, 5(2), 53. https://doi.org/10.3390/data5020053

Article Menu