Alternative cancer therapy through modeling pteridines photosensitizer quantum yield singlet oxygen production using swarm-based support vector regression and extreme learning machine

Abstract

Photodynamic cancer therapy circumvents the major side effects associated with the conventional cancer treatment methods, such as chemotherapy, surgery and exposure to radiation. Experimental measurement of photosensitizer quantum yield (PQY) singlet production of oxygen through either sensitive laser spectroscopy or luminescence detection at the wavelength of 1270 nm is costly; time consuming and intensive while unreliability of chemical traps experimental approach is of serious concern. Quantitative structure–activity relationship (QSAR) computational method proposed in the literature for computing PQY of singlet oxygen production has characteristics deviation from the measured values. PQY singlet oxygen production of twenty-nine pteridines photosensitizer compounds is modeled and predicted in this present contribution using extreme learning machine (ELM) and support vector regression (SVR) with hybridized particle swarm optimization (PSO) method for ensuring combinatory parameter selection. The performances of the developed SVR-PSO computational method are assessed using mean absolute error (MAE), correlation coefficient (CC), root mean square error (RMSE) and mean absolute percentage deviation (MAPD). The developed SVR-PSO model outperforms QSAR (2016) model with performance superiority of 34.78%, 3.65%, 17.64% and 42.16% on the basis of RMSE, CC, MAE and MAPD performance measuring parameters, respectively. The developed ELM-SINE (with sine activation function) and ELM-SIG (with sigmoid activation function) respectively outperform the existing QSAR (2016) model with improvement of 6.54% and 4.70% using R-squared metric. The demonstrated outstanding performance of the present predictive models is immensely meritorious in strengthening the potentials of alternative cancer therapy and circumventing the experimental challenges of PQY singlet oxygen production determination.

Keywords:

• Photosensitizer
• support vector regression
• pteridines
• singlet oxygen
• particle swarm optimization
• extreme learning machine

REVIEWING EDITOR:

• Jones Ian Philip PhD, Senior Editor, Imam Abdulrahman Bin Faisal University, Computer Science, Saudi Arabia
• The University of Birmingham, Metallurgy & Materials, Elms Rd, Birmingham, UNITED KINGDOM. N/A B15 2TT

SUBJECTS:

• Science
• Physical Sciences
• Material Science
• Technology
• Computer Science
• Artificial Intelligence
• Machine Learning – Design

1. Introduction

Chemical compounds that are subjected to excitation due to the incident light of a certain wavelength are referred to as photosensitizers (Frimayanti et al., 2011). Among the interesting features of photosensitizer that have profound usefulness in several applications is the singlet oxygen generation capacity (Buglak et al., 2016). Potential areas of singlet oxygen application include fine chemical synthesis, waste water treatment, solar cell and photodynamic cancer therapy among others (Buglak et al., 2016). Photodynamic cancer therapy offers a unique and novel method of cancer treatment with relatively known side effect as compared to the conventional methods, such as chemotherapy, surgery and radiation (Mfouo-Tynga et al., 2021; Wang, Liang, et al., 2021; Wang, Li, et al., 2021; Lin et al., 2021). This method of cancer treatment targets tumor cell through a combined effect of photosensitizer (which is a toxic drug that can be activated by light) and light. Photons are absorbed by a photosensitive drug and subsequently lead to drug excitation. The energy associated with the excited drug is transferred to oxygen, forming singlet oxygen otherwise known as free radicals which later oxidizes cellular structures (Oliveira et al., 2011). The transfer of energy of the excited photosensitizer is resulted from the interaction of triplet and singlet states of the sensitizer. The cancerous cell dies when the oxidative damage resulted from singlet oxygen exceeds the threshold level (Shi et al., 2015). The photosensitizer quantum yield (PQY) of singlet oxygen production is computed from the ratio of the amount of generated singlet oxygen by photosensitizer to the amount of absorbed photons. The singlet oxygen generation is influenced by the solvent polarity, oxidation potential, aggregation, energy of the triplet state, electronic configuration, steric, excimer formation and structural effects (Buglak et al., 2016; Derosa & Crutchley, 2002). Experimental determination of PQY of singlet oxygen production for photosensitizers is carried out through either sensitive laser spectroscopy or through singlet oxygen luminescence detection at the wavelength of 1270 nm. Sensitive laser spectroscopy approach is expensive, time consuming and intensive. Another Indirect experimental approach is to use singlet oxygen chemical traps. Interaction of free radicals, photosensitizer molecules and peroxides with the chemical traps limits the reliability of chemical traps approach of PQY singlet oxygen production determination while the need for the heavy water D₂O for achieving luminescence of singlet oxygen in water for luminescence detection method makes the method costly and laborious (Buglak et al., 2016). Quantitative structure–activity relationship (QSAR) has been proposed in the literature to circumvent the experimental challenges associated with PQY of singlet oxygen production while the result of the existing QSAR approach is characterized with high degree of deviation. The proposed extreme learning machine (ELM) algorithm and hybridization of support vector regression (SVR) with particle swarm optimization (PSO) have demonstrated superior performance over the existing model.

Among the well-known photosensitizers for Photodynamic cancer therapy are porfimer sodium and pteridines photosensitizer compounds. Pteridines photosensitizer is a conjugated pyrimidine and pyrazine azoheterocycles structure found in biological systems in a reduced or oxidized form Dantola et al. (2017). The pteridine compounds include the flavin (which is found in dihydro-reduced or oxidized form) and pterin sensitizer which is found in tetrahydro-reduced, half reduced or oxidized form. High potentials for electron and energy transfer, phosphorescence high yields, long lifetime excited states (triplet), photochemical activity and fluorescence characterize the oxidized forms of flavins and pterins photosensitizer. Pterins belong to heterocyclic family of compounds with six-substituted common derivatives. Pterins are classified into unconjugated and conjugated pterins in accordance to substituents functional group and molecular weight (Buglak et al., 2016). The conjugated form is characterized with short hydrocarbon chain substituents or those with one carbon atom while the unconjugated form containing larger constituents with aminobenzoic acid part. The potential of Flavins for photodynamic cancer therapy and in generating singlet oxygen is enormous and demonstrated high effectiveness as compared with the commonly used exogenous porphyrins in eradicating unwanted cells. The compositional potentials of pteridines coupled with their availability in living systems call for investigating and enhancement of their potentials in photodynamic cancer therapy using ELM and hybrid SVR and particle warm optimization algorithm presented in this work.

SVR is an intelligence computational method of solving complex real-life problems based on convex quadratic programming (Basak et al., 2007; Science et al., 2021; Drucker et al., 1996). Regularization and kernel method guide the estimation principles of this algorithm. This algorithm decreases the generalization error with the structural risk error minimization principle. This feature coupled with convergence to global solution has rendered the algorithm indispensable in many applications which include transition temperature prediction in Fe-based superconducting system (Akomolafe et al., 2021), useful life estimation of lithium-ion battery (Li et al., 2021) and energy gap modeling for enhanced photocatalytic application (Owolabi et al., 2021). The parameters in SVR algorithm can be carefully selected using different kind of approaches ranging from manual, grid-search and heuristic methods through evolutionary population-based algorithms (Wu & Wang, 2020). The evolutionary algorithm implemented in this work is PSO algorithm because of the ease of its implementation, fast convergence, high accuracy and avoidance of premature convergence (Wang, Liang, et al., 2021; Wang, Li, et al., 2021; Aboelkassem & Savic, 2021). The algorithm shares similar properties with genetic algorithm with exception of mutation and crossover operation. Hybridization of this optimization algorithm with SVR algorithm results into robust and precise model for determining pteridines PQY singlet production of oxygen. Another intelligent algorithm employed for modeling PQY singlet production of oxygen in this work is ELM (Huang et al., 2006). ELM is a single hidden layer network with intrinsic ability to randomly determine the network bias and input weights while it computes the output weights using Moore Penrose inverse approach. The uniqueness of ELM and hybrid PSO-SVR algorithms is harnessed in this work for PQY singlet production of oxygen determination.

The remaining part of this article is organized into different sections; Section 2 contains the mathematical background of the employed intelligent algorithms coupled with principles governing the implemented optimization technique. Section 3 discusses the acquisition of data for modeling and simulation. Section 4 discusses the results and compares the outcomes of the present and existing models using different performance metrics. Section 5 presents the conclusion of the manuscript.

2. Formulations of the implemented algorithms

The mathematical description of the utilized and implemented computational algorithms is discussed. This section also presents mathematical background and description of the implemented optimization algorithm.

2.1. Mathematical background of support vector regression algorithm

SVR is a learning algorithm which aims at minimizing the distance between the established fitting equation and the processed data within certain defined error limitation that is controlled by insensitive loss zone (Balogun et al., 2021; Dodangeh et al., 2020). The non-linear input vectors are transformed to higher dimensional space using kernel function while the construction of hyper-plane takes place in the feature space and controlled by the lambda hyper-parameter. Consider a given training set of data $D=\{({\eta }_1,{\varphi }_{\Delta 1}),({\eta }_2,{\varphi }_{\Delta 2}),\dots ({\eta }_m,{\varphi }_{\Delta m})\}$ where $\eta \in {\mathfrak{R}}^m,$ ${\varphi }_{\Delta }\in \mathfrak{R}$ and $j=1,2,\dots .m,$ in which ${\eta }_j$ represents the descriptors that include molecular solubility, highest occupied molecular orbital energy, dipole density, electrostatic, electronegativity and dipole while ${\varphi }_{\Delta j}$ stands for the experimentally measured quantum yield singlet oxygen production of pteridines photosensitizer. The regression equation within SVR description is shown in Equation (1)(1) $${\varphi }_{\Delta }^p={\omega }^T\eta +b$$(1) . (1) $${\varphi }_{\Delta }^p={\omega }^T\eta +b$$(1) where $\omega$ the weight vector coefficients while $b$ is the bias part

SVR learning algorithm finds a function through precise approximation, in which the expression ${\varphi }_{\Delta j}^p\approx {\varphi }_{\Delta j}$ is satisfied, where ${\varphi }_{\Delta }^p$ stands for the predicted quantum yield singlet oxygen production for pteridines photosensitizer. The algorithm aims at determining $\omega$ and $b$ through solving the optimization problem which is governed by $\epsilon$-insensitive loss function that is capable to constrain the acquired expression ${\varphi }_{\Delta }^p$ to error limit of $\epsilon .$ This objective is achieved through minimization of Euclidean norm $\frac{{\Vert \omega \Vert }^2}{2}$ subjected to the constraints presented in Equation (2)(2) $$\{\begin{array}{c}{\varphi }_{\Delta j}-({\omega }^T{\eta }_j+b)\le \epsilon \\ ({\omega }^T{\eta }_j+b)-{\varphi }_{\Delta j}\le \epsilon \end{array}$$(2) (Owolabi, 2019a, 2019b; Rui et al., 2019; Murillo-Escobar et al., 2019; Ju et al., 2019). (2) $$\{\begin{array}{c}{\varphi }_{\Delta j}-({\omega }^T{\eta }_j+b)\le \epsilon \\ ({\omega }^T{\eta }_j+b)-{\varphi }_{\Delta j}\le \epsilon \end{array}$$(2)

In a real sense, every data-points $(\eta ,{\varphi }_{\Delta })$ in $D$ cannot be approximated within $\epsilon$ error limit and can alter the feasibility of the optimization problem. This is addressed through incorporation of positive non-zero variables called slack variables $({\xi }_j,{\xi }_j^{*}).$ With this incorporation, modification of the convex optimization problem is contained in Equation (3)(3) $$\text{Minimize}\frac{{\Vert \omega \Vert }^2}{2}+C\sum _{j=1}^m\left({\xi }_j+{\xi }_j^{*}\right)$$(3) and the associated constraints are mathematically shown in Equation (4)(4) $$\{\begin{array}{l}{\varphi }_{\Delta j}-({\omega }^T{\eta }_j+b)\le \epsilon +{\xi }_j\\ ({\omega }^T{\eta }_j+b)-{\varphi }_{\Delta j}\le \epsilon +{\xi }_j^{*},\text{where}{\xi }_j,{\xi }_j^{*}\ge 0\end{array}$$(4) (Olubi et al., 2021). (3) $$\text{Minimize}\frac{{\Vert \omega \Vert }^2}{2}+C\sum _{j=1}^m\left({\xi }_j+{\xi }_j^{*}\right)$$(3) (4) $$\{\begin{array}{l}{\varphi }_{\Delta j}-({\omega }^T{\eta }_j+b)\le \epsilon +{\xi }_j\\ ({\omega }^T{\eta }_j+b)-{\varphi }_{\Delta j}\le \epsilon +{\xi }_j^{*},\text{where}{\xi }_j,{\xi }_j^{*}\ge 0\end{array}$$(4)

The parameter C in Equation (3)(3) $$\text{Minimize}\frac{{\Vert \omega \Vert }^2}{2}+C\sum _{j=1}^m\left({\xi }_j+{\xi }_j^{*}\right)$$(3) which must be greater than zero is called the regularization factor and functions in controlling the trade-off between the maximum allowable deviation and the complexity of the model. When the dual formulation is invoked with the implementation of Langrage multipliers (${\psi }_j$ ${\psi }_j^{*}$) which describe the algorithm optimization problem, the chosen kernel function $\chi ({\eta }_j,{\eta }_i)$ helps in returning pairwise dot product in high feature space without explicit data mapping. This results into expression shown in Equation (5)(5) $$\mathrm{\text{max}}-\frac{1}{2}\sum _{j,i=1}^m\left({\psi }_j-{\psi }_j^{*}\right)\left({\psi }_i-{\psi }_i^{*}\right)\chi ({\eta }_j,{\eta }_i)-\epsilon \sum _{j=1}^m\left({\psi }_j+{\psi }_j^{*}\right)+\sum _{j=1}^m{\varphi }_{\Delta j}\left({\psi }_j-{\psi }_j^{*}\right)$$(5) subject to the conditions depicted by Equation (6)(6) $$\{\begin{array}{c}\sum _{j=1}^m\left({\psi }_j-{\psi }_j^{*}\right)=0\\ {\psi }_j,{\psi }_j^{*}\in [0,C]\end{array}$$(6) (Owolabi, 2023) (5) $$\mathrm{\text{max}}-\frac{1}{2}\sum _{j,i=1}^m\left({\psi }_j-{\psi }_j^{*}\right)\left({\psi }_i-{\psi }_i^{*}\right)\chi ({\eta }_j,{\eta }_i)-\epsilon \sum _{j=1}^m\left({\psi }_j+{\psi }_j^{*}\right)+\sum _{j=1}^m{\varphi }_{\Delta j}\left({\psi }_j-{\psi }_j^{*}\right)$$(5) (6) $$\{\begin{array}{c}\sum _{j=1}^m\left({\psi }_j-{\psi }_j^{*}\right)=0\\ {\psi }_j,{\psi }_j^{*}\in [0,C]\end{array}$$(6)

The SVR approximated function in terms of the kernel and Lagrange multipliers is presented in Equation (7)(7) $${\varphi }_{\Delta }^p=\sum _{j=1}^m\left({\psi }_j-{\psi }_j^{*}\right)\chi ({\eta }_j,{\eta }_i)+b$$(7) while best function that transforms the data for PQY singlet oxygen production prediction is presented in Equation (8)(8) $$\chi ({\eta }_j,{\eta }_i)=\hspace{0.17em}\mathrm{\text{exp}}\hspace{0.17em}\left(\alpha {\left|{\eta }_j-\eta \right|}^2\right)$$(8) . (7) $${\varphi }_{\Delta }^p=\sum _{j=1}^m\left({\psi }_j-{\psi }_j^{*}\right)\chi ({\eta }_j,{\eta }_i)+b$$(7) (8) $$\chi ({\eta }_j,{\eta }_i)=\hspace{0.17em}\mathrm{\text{exp}}\hspace{0.17em}\left(\alpha {\left|{\eta }_j-\eta \right|}^2\right)$$(8) where $\alpha$ represents the parameter of the kernel function also called kernel option.

The regularization factors C, kernel option $\alpha$ and epsilon $\epsilon$ are the SVR hyper-parameters that significantly influence the model performance. These parameters are well selected and evolutionarily optimized through optimization algorithm based on swarm movement.

2.2. Particle swarm optimization algorithm

PSO technique is a meta-heuristic population-based technique of searching for global solution with unique features of simplicity, convenience and fast convergence speed (Wang, Liang, et al., 2021; Wang, Li, et al., 2021; Zou, 2021; Beheshti, 2020). The operational working principle of the algorithm involves initialization of a definite population number of random solutions where each solution (referred as particle) is allocated a velocity with which the particle searches within N-dimensional search space. The present position of the particle is controlled using particle fitness while the particle always keeps the best position in its memory (Yildiz, 2019; Yıldız et al., 2019). Consider a swarm of particles occupying N-dimensional search space where $y_k(j)=[y_{k1},y_{k2},\dots ,y_{\mathit{\text{kN}}}]$ represents the spatial position of $\mathit{\text{kth}}$ particle within the search space. Supposing $V_k(j)=[v_{k1},v_{k2},\dots v_{\mathit{\text{kN}}}]$ represents the velocity of the particles within the swarm while the best individual particle position is expressed as $r_k(j)=[r_{k1},r_{k2},\dots r_{\mathit{\text{kN}}}]$ and the best position of each of the particle is represented as $z_k(j)=[z_{k1},z_{k2},\dots z_{\mathit{\text{kN}}}].$ Equations (9)(9) $$v_k(j)=\lambda v_k(j-1)+c_{1}r\mathit{\text{and}}[r_k(j-1)-y_k(j-1)]+c_{2}r\mathit{\text{and}}[z_k(j-1)-y_k(j-1)]$$(9) and (10)(10) $$y_k(j)=v_k(j)+y_k(j-1)$$(10) , respectively, update the velocity and the position of $\mathit{\text{kth}}$ particle as the particle navigates between $j-1$ and $j$ generation (9) $$v_k(j)=\lambda v_k(j-1)+c_{1}r\mathit{\text{and}}[r_k(j-1)-y_k(j-1)]+c_{2}r\mathit{\text{and}}[z_k(j-1)-y_k(j-1)]$$(9) (10) $$y_k(j)=v_k(j)+y_k(j-1)$$(10) where $\lambda$ represents the inertial weight, $\mathit{\text{rand}}$ is the random number that spans within $[0,1]$ interval which enhances the randomness of the search, $c_1$ and $c_2$ are the inertial factors which are non-negative acceleration parameters that influences the global and local learning rate (Owolabi, 2019a, 2019b). The particle number satisfies $k=1,2,\dots .,p,$ where p is the maximum number of particles in the swarm.

2.3. Extreme learning machine

ELM is a class of neural network machine learning algorithm with single-hidden layer feedforward networks (SLFNs) (Huang et al., 2006; Oyeneyin, 2021; Alqahtani et al., 2022). Notably, error optimization using gradient decent during back-propagation (BP) is quite a time and computationally intensive process. This usually causes significantly slow training time. Also, the models generated thereof are usually non-generalized as the initial values of defined weights are context and use-case specific (Oyeneyin, 2021). In ELM, the weights between the input and the hidden layer are assigned random values, typically by using a uniform probability distribution. They remain unchanged throughout the training process, thus eliminating the need for the iterative hidden layer adjustment usually associated with neural networks. The weights between the hidden and output layer can then be derived analytically using the Moore–Penrose pseudo-inverse. The basic structure of an ELM is given in Equation (11)(11) $$\sum _{\mathrm{i}}^h{\mathrm{ɸ}}_{i\mathrm{ }}g\left(w_i.x_j+b_i\right)=y_j$$(11) . (11) $$\sum _{\mathrm{i}}^h{\mathrm{ɸ}}_{i\mathrm{ }}g\left(w_i.x_j+b_i\right)=y_j$$(11) where ${ɸ}_{i }$ is the output weight vector, $w_i$ is the input weight vector, $x_j$ and $y_j$ denote the ith input and output vector, $b_i$ is the bias of the ith hidden neuron and $g$ is the activation function of the ith neuron. Equation (11)(11) $$\sum _{\mathrm{i}}^h{\mathrm{ɸ}}_{i\mathrm{ }}g\left(w_i.x_j+b_i\right)=y_j$$(11) can be further expressed as a system of matrices as shown in Equation (12)(12) $$H.\beta =T$$(12) . (12) $$H.\beta =T$$(12)

The matrix of the output weight $\beta$ can be found analytically by solving the system of matrices in Equation (12)(12) $$H.\beta =T$$(12) . ELM has demonstrated better performance over similar neural networks with its extremely fast training time and better generalization ability (Shamsah & Owolabi, 2020).

3. Computational methodology of the hybrid model

The data acquisition strategies and description of descriptors are presented in this section. Method of algorithm hybridization is also presented.

3.1. Dataset description and acquisition

PQY singlet oxygen production for twenty-nine pteridines photosensitizer was modeled using molecular descriptors which include molecular solubility, highest occupied molecular orbital energy, dipole density, electrostatic, electronegativity and dipole. The measured dataset used for the modeling and simulation was acquired and extracted from the literature (Buglak et al., 2016). Each of the considered descriptors has direct effect on PQY singlet oxygen production for pteridines photosensitizer. The molecular solubility is inversely related with PQY singlet oxygen production. This relationship can be well understood from the deactivation induced by charge-transfer that increases with polar solvent. Thereby, the charge-transfer exciplex gains no access to water molecules due to low solubility molecules characterizing the larger part of the substituents. This consequently lowers the rate at which these compounds quench the singlet oxygen. The interaction of singlet oxygen quenching and radiation-less energy transfer due to diplole-dipole explains the influence of diplo energy on PQY singlet oxygen production. The electron attraction tendency of molecule is measured by electronegativity. In a case of the quencher molecule becoming the electron donor in charge transfer interaction, the total rate constant of singlet oxygen quenching increases which subsequently decreases the electronegativity and increases the highest occupied molecular orbital energy. Preliminary outcomes of analysis (statistical) performed on the dataset of the investigated 29 pteridines photosensitizer are presented in Table 1.

Table 1. Statistical results of the investigated pteridine photosensitizer.

	Correlation coefficient	Maximum	Mean	Minimum	Standard deviation
Highest occupied molecular orbital	−0.8023	−8.25	−9.06	−9.98	0.45
Dipole	0.1958	8.35	5.41	2.45	1.81
Molecular solubility	0.4062	−0.09	−1.72	−3.97	1.12
Electrostatic	0.4573	−0.57	−0.68	−0.76	0.04
Electronegativity	0.8588	5.65	5.10	4.41	0.37
Diploe density	0.4960	0.16	0.08	0.02	0.04
Quantum	1.0000	1.73	1.17	0.30	0.50

The dipole, electrostatic and dipole density are weakly correlated with the quantum yield of singlet oxygen production while the highest occupied molecular orbital shows negative correlation. The coefficient of correlation shows the extent and degree of linear connection between the descriptors and targets. The presented maximum and minimum values presented in table shows the range of the dataset contents while the consistency in measurement is measured by the presented standard deviations. The uniqueness and the relevancy of the chosen descriptors enhance the developed predictive models for precise measurement of quantum yield of singlet oxygen production.

3.2. Computational hybridization of PSO with SVR algorithm

Modeling of quantum yield of singlet oxygen of pteridines photosensitizer for photodynamic cancer therapy application was carried out using computing part of MATLAB. Randomization of the available descriptors (molecular solubility, highest occupied molecular orbital energy, dipole density, electrostatic, electronegativity and dipole) and their corresponding targets was initiated prior to dataset partitioning. The essence of randomization is to ensure even and uniform data-points distribution so that the developed model is not biased toward the training or testing phase. The twenty-nine pteridines compounds available for simulation were separated into training and testing set in the ratio of 8:2. The training set of data was employed for support vector acquisition while the testing data samples were utilized for validating the performance of the acquired support vectors during training phase. The step-by-step procedures of algorithm hybridization are itemized as follows:

• Step I: Particle randomization and initialization: at $j=0,$k-number of particles were randomized uniformly and initialized $r_{\mathit{\text{kN}}}(0)$ within $[r_{\mathrm{\text{max}}},r_{\mathrm{\text{min}}}]$ such that $r_{\mathrm{\text{min}}}=-r_{\mathrm{\text{max}}}.$ Initialization of the velocity take similar approach with uniform distribution of particle characterized with $v_{\mathit{\text{kN}}}(0)$ within $[v_{\mathrm{\text{max}}},v_{\mathrm{\text{min}}}]$ range in which $v_{\mathrm{\text{min}}}=-v_{\mathrm{\text{max}}}.$ The maximum velocity in the given direction is expressed as $v_k{}_{\mathrm{\text{max}}}=\frac{y_k{}_{\mathrm{\text{max}}}-y_k{}_{\mathrm{\text{min}}}}{i}$ with $i$ number of the interval chosen.

• Step II: Fitness determination through objective function implementation: Each of the particle within the swarm was set at its current best position $z_k^{*}(0)=z_k(0)$ while the corresponding fitness function is represented as $F_k^{*}(0)=F_k(0).$ The particle with best fitness $F_{\mathit{\text{best}}}$ (corresponding to lowest root mean square error [RMSE]) was set as global fitness with position $z_k^{**}(0)$ while $F_k^{**}(0)=F_{\mathit{\text{best}}}.$ The following steps were implemented while determining the fitness of each of the particle in the swarm.

• Step A: Mapping function: A function was selected among several probable functions that can transform the input set of data to high dimensional space. The available functions include linear, sigmoid, Gaussian and polynomial.

• Step B: Selection of a particle in a swarm for fitness computation: a particle among the swarm which houses the probable values of the hyper-parameter in the order of penalty factor, kernel parameter and epsilon, was selected. The selected particle together with the training set of data and the chosen kernel function were utilized for SVR training and the trained model was evaluated using RMSE-training. The acquired support vectors were saved for future use.

• Step C: Validation of the fitness of the trained model: Using the testing data samples alone alongside with the acquired support vectors in Step B, the fitness of the trained model was validated using RMSE-testing.

• Step D: Repeat Step B and Step C for the entire particle in the swarm and save the support vectors together with the values of the parameters that measure the fitness of each of the particle.

• Step III: Iteration update: the subsequent iteration was updated through $j=j+1$ implementation.

• Step IV: Updating inertial weight $\lambda$: The inertial weight is updated in accordance to $\lambda (j)=\beta \lambda (j-1)$ where $\beta$ represents a constant which approaches unity in value.

• Step V: Update in velocity and position: Velocity and position of each of the particle are updated regularly in every iteration until stopping condition is attained using Equations (9)(9) $$v_k(j)=\lambda v_k(j-1)+c_{1}r\mathit{\text{and}}[r_k(j-1)-y_k(j-1)]+c_{2}r\mathit{\text{and}}[z_k(j-1)-y_k(j-1)]$$(9) and (10)(10) $$y_k(j)=v_k(j)+y_k(j-1)$$(10) , respectively.

• Step VI: Updating individual best position if necessary: Supposing $F_k(j)<F_k^{*}(j),$ the individual best position as well as the corresponding fitness is updated as $z_k^{*}(j)=z_k(j)$ and $F_k^{*}(j)=F_k(j).$ If $F_k(j)<F_k^{*}(j)$ condition is not met, proceed without updating.

• Step VII: Update the global best position if necessary: if $F_{\mathrm{\text{min}}}(j)<F^{**}(j)$ condition is met, then $F^{**}(j)=F_{\mathrm{\text{min}}}$ and $z^{**}(j)=z_{\mathrm{\text{min}}}(j).$

• Step VIII: Stopping condition: Check if $F_{\mathrm{\text{min}}}=0$ or $F_{\mathrm{\text{min}}}$ attains the same value for fifty consecutive iterations. Otherwise, return to Step III.

4. Results and discussion

The outcome of the simulation carried out in this research work is presented. The results of the optimization as well the performance comparison between the present and existing models are also presented.

4.1. Convergence of the hybrid model at different swarm population

The variation of the convergence of the developed hybrid SVR-PSO model to the number of swarm exploring and exploiting the search space is shown in Figure 1. For the error convergence presented in Figure 1(a), different number of swarm population has different starting point of convergence, however, they all show similar global solution after several iterations. It indicates high level of robustness of the SVR-PSO model. The convergence of penalty/regularization hyper-parameter factor at various number of swarm population is presented in Figure 1(b). The exploitation and exploration ability of the swarm within the defined search space are well investigated. The convergence of the regularization factor when 50 swarms explored the search space appears faster than the case when 100 and 500 swarms were exploiting and exploring the space. Figure 1(c) presents the convergence of kernel parameter for different size of swarm in the search space at various number of iteration. The model converges faster when the number of swarm in the population was set at 50 as compared to the cases when the number of swarm was maintained at 100 and 500. Similar convergence pattern is obtained for the variation of epsilon at different iteration stages as presented in Figure 1(d).

Figure 1. Convergence of the developed SVR-PSO model at various iteration for different number of swarms (a) root mean square error; (b) regularization factor; (c) kernel parameter; (d) epsilon.

The obtained values of the optimum hyper-parameters as extracted from the PSO are presented in Table 2.

Table 2. Results of the optimization of SVR hyper-parameters using PSO.

Hyper-parameters	Optimum value
Regularization factor	167.7855
Epsilon	0.1934
Hyper-parameter lambda	1.00E-07
Kernel parameter	0.6303
Number of swarm in population	50
Kernel function	Gaussian

4.2. Comparison between developed intelligent models

The performance of the developed three models is compared on the basis of R-squared at different stages of model development and presented in Figure 2. During pattern acquisition, the developed PSO-SVR model demonstrated superior performance over the other two models (ELM-SINE and ELM-SIG) followed by ELM-SINE and ELM-SIG model. While implementing the acquired patterns for validation set of data, ELM-SIG model shows superior performance followed by ELM-SINE and PSO-SVR model. Table 3 presents the values of each of the assessment parameters for different stages of model development. The assessment parameters include the correlation coefficient (CC), root mean square error (RMSE), mean absolute error (MAE) and R-squared.

Figure 2. Performance comparison between developed models using R-squared parameter (a) training stage and (b) testing stage.

Table 3. Performance assessment parameters for the developed models.

	Training set				Testing set
	CC	R-squared	RMSE	MAE	CC	R-squared	RMSE	MAE
ELM-SINE	0.9874	0.9750	0.1498	0.1107	0.9894	0.9790	0.1803	0.1278
ELM-SIG	0.9742	0.9491	0.2134	0.1559	0.9988	0.9976	0.0913	0.0758
PSO-SVR	0.9946	0.9892	0.1795	0.1717	0.8601	0.7398	0.1910	0.1706

The weights associated with the developed ELM-based models are presented in Table 4. The model weights contain the randomly generated biases, computed output weights using Moore–Penrose pseudo-inverse approach and optimized input weights corresponding to input descriptors.

Table 4. Weights associated with the developed ELM-based models.

ELM-Sine	j	bias	Output weights	First input weight	Second input weight	Third input weight	Forth input weight	Fifth input weight	Sixth input weight
	1	−0.77543	0.537917	−0.15378	0.785845	−0.85201	0.19498	0.935731	−0.98996
	2	0.568856	−0.20799	0.311146	0.406446	0.368192	−0.32938	0.457886	1.730396
	3	−0.41686	−0.45412	0.445845	0.111476	−0.19522	−0.40155	0.240478	0.511089
	4	0.207067	−0.92553	0.062419	−0.63113	0.96567	−0.09481	0.763898	−1.53541
	5	0.928845	0.34659	−0.78236	−0.57594	−0.19563	−0.15471	0.759327	0.575405
	6	−0.13503	−0.14087	0.263533	−0.84531	0.241344	−0.28079	0.740648	1.291197
	7	0.389504	−0.09652	−0.747	0.827601	−0.69126	0.116638	0.743688	−0.21995
	8	0.516199	0.219714	−0.73139	0.41343	−0.23731	0.485091	0.10592	−2.51514
	9	−0.13472	−0.88119	−0.80281	0.115578	−0.67773	−0.15133	0.68156	−1.11964
	10	0.310996	−0.36838	−0.71595	−0.37314	0.516225	−0.14129	0.463261	−0.58486
	11	−0.78049	0.545444	−0.6635	−0.66759	0.742222	−0.75025	0.212163	0.057294
	12	0.86752	0.392866	−0.6075	0.244995	−0.29845	−0.95113	0.098519	−0.23167
	13	−0.62508	−0.74934	−0.36504	0.975869	0.371071	−0.41963	0.823574	−0.33911
	14	−0.46764	−0.7397	−0.36714	−0.65914	−0.4117	−0.36496	0.17501	−1.89123
	15	0.595661	−0.8153	−0.56487	−0.48442	0.061259	0.30738	0.16357	−0.9823
	16	−0.02479	−0.98436	−0.49792	−0.2064	0.664847	0.913872	0.665987	0.866364
ELM-Sig	1	−0.0924	−0.20148	−0.78757	−0.15433	0.333056	−0.26217	0.120612	−53.6521
	2	−0.13522	0.053752	−0.25518	0.095742	−0.64374	−0.07855	0.589507	−6.62629
	3	0.650628	−0.1664	−0.60376	0.885474	−0.74397	0.963276	0.226188	−5.96943
	4	−0.83306	0.31372	−0.02062	−0.16451	0.998161	−0.68719	0.384619	559.6129
	5	−0.73366	0.255947	−0.32101	0.966105	−0.65776	0.711046	0.582986	−1207.89
	6	−0.65322	−0.41603	0.903261	−0.39709	−0.9348	0.289529	0.251806	101.2744
	7	−0.21812	−0.1367	0.840664	0.402198	0.1224	−0.24746	0.290441	3413.36
	8	0.662759	−0.96903	−0.89465	0.332678	0.763733	−0.61815	0.617091	48.17706
	9	0.606729	0.968127	0.475716	0.078253	0.338351	−0.14349	0.265281	−16.6354
	10	−0.87906	−0.66566	−0.46176	0.396211	−0.61913	−0.03596	0.824376	−2857.74

4.3. Performance of the developed predictive models and comparison with the existing QSAR (2016) model

The performances of the developed hybrid SVR-PSO and ELM-based models are evaluated using CC, R-squared, RMSE and the MAE and presented in Table 5. The values of CC, R-squared, RMSE and MAE associated with the developed SVR-PSO model are 0.9908, 0.9816, 0.1815 and 0.1715, respectively, while the computed assessment parameters for the developed ELM-SIG model are 0.9781, 0.9566,0.1978 and 0.1421, respectively. For the developed ELM-SINE model, the values of CC, R-squared, RMSE and MAE obtained are 0.9866, 0.9735, 0.1555 and 0.1137, respectively. Each of the developed models has higher values of CC and R-squared and lower values of errors (RMSE and MAE) as compared with the existing model. On the basis of comparison of the developed SVR-PSO with the existing QSAR (2016) (Buglak et al., 2016) model, performance improvement of 34.78%, 17.64%, 3.65% and 7.44% were, respectively, obtained on the basis of RMSE, MAE, CC and R-squared performance metrics.

Table 5. Comparison between the models and the performance improvement of SVR-PSO over QSAR (2016) (Buglak et al., 2016).

Performance assessment parameter	QSAR (2016) (Buglak et al., 2016)	SVR-PSO (this work)	ELM-SIG (this work)	ELM-SINE (this work)	% improvement of SVR-PSO over QSAR (2016)	% improvement of ELM-SIG over QSAR (2016)	% improvement of ELM-SINE over QSAR (2016)
MAE	0.2083	0.1715	0.1421	0.1137	17.6440	31.7832	45.4304
RMSE	0.2783	0.1815	0.1978	0.1555	34.7797	28.9255	44.1222
CC	0.9558	0.9908	0.9781	0.9866	3.6547	2.3247	3.2227
R-squared	0.9136	0.9816	0.9566	0.9735	7.4431	4.7034	6.5493

Comparison of the developed ELM-SIG model with the existing QSAR (2016) model (Buglak et al., 2016) shows the superiority of the developed ELM-SIG model with enhancement of 31.78%, 28.93%, 2.32% and 4.70% using MAE,RMSE, CC and R-squared performance assessment, respectively. Similarly, the developed ELM-SINE model outperforms the existing QSAR (2016) model (Buglak et al., 2016) model with improvement of 45.43%, 44.12%, 3.22% and 6.54% using performance metrics of MAE, RMSE, CC and R-squared, respectively. Cross-plot correlation between the estimated and measured quantum yield of singlet oxygen production is presented in Figure 3 for the present and existing models.

Figure 3. Correlation cross-plot between the measured and estimated quantum yield of singlet oxygen production.

The data-points for the investigated pteridine photosensitizer estimated by the developed SVR-PSO model, ELM-SIG model and ELM-SINE model are well aligned with the measured values while the deviation in the estimates of QSAR (2016) model (Buglak et al., 2016) can be inferred from its non-aligned data-points.

5. Conclusion

The quantum yield of singlet oxygen production for 29 different compounds of pteridine photosensitizer is modeled using ELM and hybridization of SVR with PSO algorithm. The predictors to the developed model include molecular solubility, highest occupied molecular orbital energy, dipole density, electrostatic, electronegativity and dipole while the developed model was validated using experimentally measured quantum yield of singlet oxygen production. The performance of the developed SVR-PSO model is compared with the existing QSAR (2016) model in the literature on the basis of RMSE, CC, MAE and MAPE and the hybrid model shows superior performance over the existing model. Comparison of the developed ELM-SIG model with the existing QSAR (2016) model shows the superiority of the developed ELM-SIG model with enhancement of 31.78%, 28.93%, 2.32% and 4.70% using MAE, RMSE, CC and R-squared performance assessment, respectively. Similarly, the developed ELM-SINE model outperforms the existing QSAR (2016) model with improvement of 45.43%, 44.12%, 3.22% and 6.54% using performance metrics of MAE, RMSE, CC and R-squared, respectively. The accuracy, robustness and outstanding performance of the developed predictive models would definitely promote quick, easy and less costly method of determining quantum yield of singlet oxygen production for photodynamic cancer therapy treatment. The developed model is limited to pteridine PQY of singlet oxygen production.

Disclosure statement

Author has no competing interest.

Data availability statement

The raw data needed to reproduce these findings are cited in Section 3.1 of the manuscript.

Additional information

Notes on contributors

Nahier Aldhafferi

Dr. Nahier Aldhafferi is an associate professor at the Department of Compteure Information Systems, College of Computer Science and Information Technology, Imam Abdurrahman Bin Faisal University in Saudi Arabia. He received his BS in 2005 from Dammam Teachers College, Saudi Arabia. He also received his master’s degree in 2009 in Internet Technology from Wollongong University, Australia. Currently, he has a PhD in Information Technology from New England University, Australia. He has published alot of research in Data Science and Information Privacy and Artificial Intelligence. His current research interests are in the fields of Data Science, Data Mining, E-gov Services Integration, Privacy, Machine Learning, and Software Engineering. He can be contaced at email: [email protected]