Stroke is the second most prevalent cause of death [1], hence the special relevance of effective pharmacophores for its treatment. Development and bench-to-bedside translation of pharmaceuticals is complex and lengthy, involving multiple preclinical and clinical steps. Optimization of such schemes by means of mathematical modeling using computer technologies is a well-defined task of applied bioinformatics.

Translational algorithms of computational pharmacology are increasingly recognized as integral tools in preclinical and clinical trials for innovative drugs. However, the expert judgments on this essential trend are non-uniform and depend on social/ cultural background and economical/legal context.

New nanoparticle-based pharmacophores are a long-term focus of preclinical studies programs passed by the Shanghai Cooperation Organisation (SCO) for 2020–2025 [2, 3].

Significantly, this regional political and economical association of nine countries (China, Russia, Iran, India, Pakistan, Kazakhstan, Uzbekistan, Tajikistan and Kyrgyzstan) with a total budget of almost 26.8 billion USD currently provides support for ongoing developments in the field of preclinical studies optimization by emphasizing the influence of nanopharmacology models in silico [3, 4].

The apparent ethical and economical discrepancies between 'western' and 'eastern' perspectives on research, development and registration of new drugs [2, 3] signal an opportunity of modifying our attitude towards in silico methods in order to reduce both the costs and the lengths of pharmacophore translation. Sure enough, financial limitations and related obstacles [4] dictate a carefully weighted position with regard to a number of sensitive issues concerning translational medicine. On the other hand, we encounter the amassing experience of preclinical research optimization in Eurasian countries over the last 5–7 years.

For instance, recent translational projects controlled by SCO encompass current research on the neuroprotective potential of porphyrin-fullerene nanoparticles releasing ²⁵Mg²⁺. A remarkable aspect of this endeavour is the active use of in silico platforms designed to analyze allometric and in vitro data to predict the numerical values and variability levels for the most relevant pharmacokinetic and pharmacodynamic parameters [5–7].

After almost 15 years of studies on pharmacological effects of paramagnetic metal isotopes [4–11], three of them, ²⁵Mg, ⁴³Ca and ⁶⁷Zn, have been qualified as promising candidates for targeted delivery with the use of porphyrin-C60-based nanocarriers [8–11]. The above-mentioned reasons necessitate the use of in silico platforms with corresponding translational scenarios as appropriate [4, 5, 11, 12].

Modified algorithms of non-Markovian population dynamics, along with computational models based on ordinary differential equation (ODE) systems, are effective for predicting selective accumulation ²⁵Mg, ⁴³Ca, ⁶⁰Co and ⁶⁷Zn divalent metal isotopes released from PMC16 porphyrin-fullerene nanoparticles in brain cells [3].

This study aimed to develop a drug-oriented mathematical model in the ordinary differential equation system format for in silico implementation in MATLAB and adjust its adequatness using in vivo rat model.

METHODS

In silico model

The ODE system (1) used as a basis for the mathematical model implemented in silico was built analytically using published biomedical evidence and our own experimental data concerning the in vivo effect of ²⁵Mg²⁺-releasing porphyrinfullerene nanoparticles (NPs) on the cerebral ischemic stroke pathogenesis. The code was written in MATLAB and implemented in MATLAB/SIMULINK version 2021b. The ODE system was solved using the Runge–Kutta midpoint method with the ode45 solver; T_1/2 = 9.0 h, T_max = 25 h, С₀ = 62 µg/mL.

The behavior of the model in the parameter space was studied based on the theory of nonlinear dynamic systems with the method known as 'parametric analysis' (or 'bifurcation analysis' in synergetics) [14]. The method is effective for studying complex nonlinear processes in non-equilibrium systems of various types: physical, chemical, social, biological, etc., using their mathematical models on a computational experiment basis [14].

The mathematical model was built in the ODE system format by applying the five-compartment framework and choosing parameters to satisfy the in vivo experimental data. After reduction to a Cauchy problem by numerical method, the model was translated into machine code and further investigated by the parametric analysis of nonequilibrium processes and systems to bring it into line with the available experimental data.

The parametric analysis method [14] was further applied for studying the behavior of the model under various dynamic combinations of parameters to test the validity and adequacy of the proposed in silico platform, its power in solving pharmacokinetic problems, as well as the prospects of its use in optimization of the algorithms for preclinical studies.

In vivo experiment

NPs

The water-soluble PMC16-RX samples were courtesy of Dr N. Amirshahi, Amirkabir University of Technology, Tehran, Iran.

Animals

Wistar Albino Glaxo male rats (n = 30, purchased from BioPitomnik STEZAR facilities; Vladimir, Russia) of 180–220 g body weight were maintained on a standard fortified diet. The animals were fasted for 24 h before the experiment. Each experimental point involved 3 animals; each measurement was performed in 5–6 replicates. All manipulations involving animals (including withdrawal from the experiment) were performed in accordance with ARRIVE guidelines.

NP administration

NPs were administered to rats at 1.0 mg/kg and/or 20.0 mg/kg as a single intravenous injection using 15 mM Tris-HCl (рН 7.80) as a solvent. The animals were decapitated 12 h post-injection and brain tissue samples were homogenized in 5–7 volumes of 20 mM Трис-HCl (рН 8.0)/10 mM MgCl₂/1.5 mM NaCl/2.0 mM EDTA/25 mM sucrose/2.0% Triton X-100 (w/w) using a glass/ teflon Potter homogenizer at 1800 rpm (+4 °C).

Brain tissue homogenate fractionation

To isolate the S125 cytosolic fraction, the carefully washed homogenates were ultracentrifugated at 125,000 g for 4 h at +4 °C in a Spinco L5-65B ultracentrifuge (Beckman; USA), rotor SW 27.1. The S125 supernatants were collected and the protein content was measured by the conventional colorimetric Bradford assay.

The S125 aliquots were mixed with 10 volumes of icecold aceton for subsequent overnight incubation at +4 °C. The precipitates were pelleted by centrifugation (20,000 rpm, 20 min, +4 °C) and discarded; the supernatants were collected and preserved for UV-spectrometry and CZE analysis.

To assess the NP extractability from biomaterial the aceton-insoluble pellets were dried at 25 °C and dissolved in 15 mM ammonium sulfate (рН 8.80)/0.1% SDS/2.5 mM EDTA/1.0% 2-mercaptoethanol (20 : 1, v/w) for subsequent 60 kHz ultrasonication at 40 °C for 60 min. The samples were subsequently analyzed by CZE (see below).

Ischemic stroke in vivo model

The stroke was modeled using conventional technique of of the middle cerebral artery filament occlusion described elsewhere [15].

Original protocol of the capillary zonal electrophoresis (CZE)

The aceton-soluble S125 extracts were concentrated in a rotary evaporator to 0.2–0.3 mL volume and diluted with 30 mM ammonium phosphate (рН 8.80) 25 : 1, v/v. A 10 µl aliquot of this sample was injected into the receiver of a P/ACE MDQ Plus CZE analytical system (ALGIMED; Belarus) equipped with a UV-VIS 770 KS online fraction detector with a 440 nm monochromatic filter (Prince Technologies; Netherlands).

The CZE process was initiated and accomplished over 10 min at +6°C in a 50 µm (diameter) × 75 mm (effective length) quartz capillary containing UV-transparent silica gel saturated with SJX40 electrolyte pH 8.0 (SCLEX, BV; Netherlands).

Fractionation mode: 115 V/60 Hz/300 W/capillary. Block of planigraphic/temporal data analysis DAX DATE 220 LK (SCIEX BV; Netherlands). The calibration was carried out conventionally, using an internal standard in accordance with the algorithm recommended by the analytical system manufacturer (ALGIMED; Belarus); the calibration data are given in tab. 1.

Pharmacokinetic distribution curves for PMC16 for different organs are given below; the confidence intervals were calculated using Student's method.

RESULTS

Our model is based on the common assumption that a nanoagent spreads throughout the body with the bloodstream, entering the heart, the liver, and the brain, and is eliminated naturally afterwards. At the same time, the post-ischemic inflammation area of the brain accumulates the nanoagent faster, which is explained by increased permeability of vascular walls in this area due to the loosening of intercellular adhesions in the endothelium. The model is represented by ODE system (1) and illustrated in fig. 1.

In ODE system (1) C_b, C_i, C_h, C_is,  and C_br are concentrations of the nanoagent in the blood, the liver, the heart, intercellular spaces of the brain and brain cells, respectively. The 'intercellular spaces' compartment was added into consideration to emphasize the striking difference in permeability of the blood brain barrier for ischemic lesions compared with non-affected brain regions. The K_is constant is adjusted with regard to the size of ischemic lesion: the bigger the inflammation focus, the higher its value.

The K_в constant reflects the nanoagent elimination. The rest of the constants reflect the rates of transition of nanoagent molecules between the compartments. It is easy to see that summation accross the system leaves us with the term K_вC_b only, i.e. the nanoagent is completely eliminated from the body. For simplicity, the removal of the nanoagent through the liver is also included in K_в, but the metabolism takes time, which is reflected by the constants corresponding to the liver.   

formula

Knowing the half-life T_1/2 [5], we can estimate the elimination constant K_в as ln2/T_1/2.

We modeled a single intravenous injection of the substance at a 0.2 mg/µg dose. Although the uniform distribution of a pharmacophore from the injection site throughout the body is known to take time, we excluded the 'diffusion' factor from the model to simplify calculations, that is, assumed the nanoagent concentration uniformity immediately after administration.

The in silico modeling curves in comparison with in vivo experimental data are shown in fig. 2.

The pharmacokinetics of the studied pharmacophore was simulated in silico over 24 h known to be the critical period of ischemic stroke development. With an increase in K_h+, the peak pharmacophore concentration in the myocardium shifts to the left and increases. According to published evidence [16], pharmacophore concentration rapidly peaks in the liver (within one hour) and declines rapidly, i.e. the liver clearance constant exceeds the absorption constant; at that, K_l+ is smaller than K_h+. Upon a decrease in the clearance constant K_l–, the peak concentration of pharmacophore increases sharply. The corresponding constants describing pharmacodynamics for the brain require additional data for proper validation; nevertheless, fig. 2 shows that pharmacophore accumulation in the brain is lower compared with the liver and the myocardium, which is consistent with both the literature [16], and our own in vivo experiments. Fig. 2 (fig. 2) shows the calculated and experimental pharmacokinetic curves for PMC16 in the brain, the liver and the heart. The calculated curves (shown in black) have a classical shape and numerically fit into the confidence intervals of the experimental data collected in vivo for rat model, which proves the consistency between the mathematical model and the natural experiment.

Thus, we demonstrate experimentally that РМС16 crosses the blood-brain barrier to be internalized by the post-ischemic brain cells, delivering ²⁵Mg²⁺ ions that stimulate ATP hypersynthesis and thereby exert strong neuroprotective effects. At the same time, quantitative content of the studied pharmacophore in the brain, the heart and the liver differs. Our in silico model retrieves the same results, which confirms its prognostic potential and prospects in optimization of preclinical studies algorithms for porphyrin-fullerene nanocationites as promising neuroprotectants.

The distribution of PMC16 NPs in brain regions with regard to ischemic injury is given in tab. 2.    

DISCUSSION

The pharmacokinetics of porphyrin-fullerene nanoagents for preclinical research optimization is solvable in silico on a computational experiment basis [17, 18]. However, for comprehensive optimization of research algorithms, pharmacokinetics (PK) (i.e. computational output for the pharmacophore delivery to a particular organ, as has been done by us for РМС16 in this study) should be complemented by pharmacodynamics (PD) (i.e. a corresponding output for the development of therapeutic effect towards the studied pathology, e.g. ischemic stroke). We therefore consider building a complex PK/PD model as the ultimate prospective goal of our further research in silico.

Meanwhile, such studies are complicated not only by the scarcity of comprehensible mathematical models of ischemic stroke suitable for this purpose, but also by the lack of pharmacokinetic models for innovative pharmacophores, which may show distinct drug-specific pharmacokinetic patterns, especially under conditions of targeted delivery. The rapidly expanding field of nanopharmacology offers quite a number of such agents with excellent neuroprotective properties, albeit with pharmacokinetic patterns distinct from those of conventional pharmaceuticals [19].

These considerations are immediately applicable to the domestically produced РМС16, range regarded extraordinary promising in terms of ischemic stroke therapy. The invention of these pharmacophores was based on the fundamental discovery of the magnetic isotope effect exerted by ²⁵Mg²⁺. The magnetic isotope hyperactivates the magnesium-dependent ATP syntheses in a cell; moreover, such energy metabolism hyperactivation requires a tiny amount of these ions and can proceed even in the absence of oxygen, under conditions of profound tissue hypoxia [16]. These properties are exemplified by PMC-16 — a pharmaceutical nanoagent comprising a porphyrincontaining fullerene 'sphere' C60 (porfyllerene-MC16) [20, 21].

The key issues in experimental validation of such in silico paradigm as ours include crossing of the blood-brain barrier by PMC16 and its capture (internalization) by brain cells. The task was accomplished (fig. 2; tab. 1 and tab. 2) to find both conditions fulfilled. Noteworthy, brain regions affected by the stroke were accessible to PMC16 as well (tab. 2).

Thus, the relatively low 'mass quantity level' of NP capture by rat brain cells (tab. 1) may have no correlation with the expected pharmacological effect of the agent, as the latter invariably involves the excessive ATP synthesis as a direct consequence of the ²⁵Mg²⁺ magnetic isotope effect phenomenon. At the same time, internalization of pharmacophore by target cells is definitely a priority in the advanced drug development.

A plain single-step ultracentrifugation (105,000–150,000 g) of mammalian tissue homogenates, usually pretreated with Triton X-100, ensures separation of cytosol from organelles, ribosomes, ribosomal subunits and membranous debris. The S125 fraction, a totality of soluble cytosolic compounds, is an excellent substrate to be used in pharmacophore internalization studies [3, 4, 10].

Although the obtained results are valuable as such, it should be noted that resolution and sensitivity of our CZE procedure are high enough to enable determination of the low-rate capture of PMC16-RX by rat brain, reaching estimated 4.0–8.0 ng per 1.0 mg of S125 total protein. Such content of NP, recorded 12 h after single intravenous injection of the drug, is distinctly above the background set by calibration (tab. 1), which clearly indicates the penetration of the studied xenobiotic through the blood-brain barrier (tab. 1 and tab. 2).

CONCLUSIONS

The developed in silico model describes the pharmacokinetics of РМС16 nanoparticles as promising neuroprotectants for the use in ischemic stroke therapy. The distribution of PMC16 in organs and tissues complies with a five-compartment mathematical model represented by an ordinary differential equation system and implemented computationally. The model can be fitted to specific properties of biological objects through the adjustment of numerical constants (parameters) of the model on the basis of in vivo experiments. The developed model is consistent with published evidence and newly collected experimental data and can be applied to optimize preclinical research scenarios for medical nanocationites РМС16 as neuroprotectants in ischemic stroke therapy.