Pauline Dupuis

Madani François

Thomas Baudier

Marie- Pierre Sunyach

Tristan Brahmi

Myriam Ayadi

Marie- Claude Biston

Evaluation of a dedicated software for semiautomated VMAT planning of spine Stereotactic Body Radiotherapy (SBRT)

Purpose

To determine whether SBRT of spinal metastasis using a dedicated treatment planning system (TPS) and delivered with a gantry-based LINAC could provide plans of similar quality to the Cyberknife technology. Additional comparison was also done with other commercial TPS used for volumetric modulated arc therapy (VMAT) planning.

Materials and Methods

Thirty Spine SBRT patients, previously treated in our institution with CyberKnife (Accuray, Sunnyvale) using Multiplan TPS, were replanned in VMAT with an dedicated TPS (Elements Spine SRS, Brainlab, Munich) and our clinical TPS (Monaco, Elekta LTD, Stockholm), using exactly the same arc geometry. The comparison was done by assessing differences in dose delivered to PTV, CTV and spinal cord, calculating modulation complexity scores (MCS) and performing quality control (QA) of the plans.

Results

Regardless of the vertebra level, in general, no statistical difference was found in PTV coverage between all TPS. Conversely, PTV and CTV D50% were found significantly higher for the dedicated TPS compared to others. In addition, the dedicated TPS also resulted in better gradient index (GI) than clinical VMAT TPS, whatever the vertebral level, and better GI than Cyberknife TPS for the thoracic level only. The D2% to the spinal cord was generally significantly lower with the dedicated TPS compared with others. No significant difference was found in the MCS between both VMAT TPS. All QA were clinically acceptable.

Conclusion

The Elements Spine SRS TPS offers very effective and user-friendly semi-automated planning tools and is secure and promising for gantry-based LINAC spinal SBRT.

Introduction

For oligometastatic patients, spine stereotactic body radiation therapy (SBRT) is an extending alternative to the classical palliative fractionation regimen [START_REF] Sprave | Randomized phase II trial evaluating pain response in patients with spinal metastases following stereotactic body radiotherapy versus threedimensional conformal radiotherapy[END_REF]. The improvements in systemic therapy lines allow a benefit in patient survival, and SBRT can lead to a durable local control with very limited neurologic side effects [START_REF] Husain | Stereotactic body radiotherapy for de novo spinal metastases: systematic review[END_REF], [START_REF] Sahgal | Spinal cord dose tolerance to stereotactic body radiation therapy[END_REF]. Spine SBRT planning is challenging because of the close vicinity between the spinal cord, having strict dose constraints, and the Planning Target Volume (PTV). Consequently, a high dose gradient is required to achieve the dosimetric objectives. For this reason, spine SBRT has historically been preferentially delivered using robotic radiosurgery. The CyberKnife (Accuray, Sunnyvale) provides steep dose gradients by using multiple non-coplanar beams, and can therefore be considered as the "gold standard" technology.

In parallel, new perspectives of treatments have emerged with last generation gantrybased linear accelerators (LINAC) equipped with highly performing imaging devices, thus enabling an increased accuracy in dose delivery. Therefore, dose planning and delivery of spine SBRT is now commonly performed with versatile treatment planning systems (TPS) using volumetric modulated arc therapy (VMAT) technique [START_REF] Nalichowski | Single fraction radiosurgery/stereotactic body radiation therapy (SBRT) for spine metastasis: A dosimetric comparison of multiple delivery platforms[END_REF], [START_REF] Moustakis | Treatment planning for spinal radiosurgery: A competitive multiplatform benchmark challenge[END_REF].

Simultaneously, specialized planning solutions are also elaborated. The Elements TPS (Brainlab, Munich) has been developed according to this philosophy, with specific modules dedicated to stereotactic radiosurgery (SRS): Cranial SRS, Multiple Brain Mets SRS or Spine SRS [START_REF] Trager | Evaluation of elements spine SRS plan quality for SRS and SBRT treatment of spine metastases[END_REF], [START_REF] Deshazer | A dosimetric comparative analysis of Brainlab elements and Eclipse RapidArc for spine SBRT treatment planning[END_REF]. Elements Spine SRS (E-SRS) allows an automated VMAT planning process, based on the use of templates defining beam geometry and dosimetric objectives on organs at risk (OAR) and target volumes. The optimization algorithm in E-SRS has been specifically designed to facilitate high dose gradient at the spinal cord to PTV interface. This is particularly important for targets with complex and concave shapes, e.g. when the whole vertebrae has to be treated. The optimization algorithm includes an arc duplication option enabling PTV splitting into several low complexity elementary volumes treated separately with one or several arc each. Such a decomposition should also result in treatment plans with the same or higher quality than non-dedicated TPS [START_REF] Trager | Evaluation of elements spine SRS plan quality for SRS and SBRT treatment of spine metastases[END_REF]. Another advantage of E-SRS is that the planning process is semi-automated. Hence, it should allow less inter-operator variability than conventional TPS while saving planning time.

The dosimetric performances of E-SRS TPS compared with other clinical TPS were evaluated in 4 different studies [START_REF] Trager | Evaluation of elements spine SRS plan quality for SRS and SBRT treatment of spine metastases[END_REF], [START_REF] Deshazer | A dosimetric comparative analysis of Brainlab elements and Eclipse RapidArc for spine SBRT treatment planning[END_REF], [START_REF] Piliero | Plan quality and consistency in spine radiosurgery treatment planning: comparison between automatic treatment planning with Elements Spine SRS and manual inverse planning with Varian Eclipse[END_REF], [START_REF] Saenz | A dosimetric analysis of a spine SBRT specific treatment planning system[END_REF]. In two recent studies, the dosimetric superiority of Elements Spine SRS compared with Eclipse TPS (Varian Medical Systems, Palo Alto, USA), in terms of spinal cord sparing [START_REF] Trager | Evaluation of elements spine SRS plan quality for SRS and SBRT treatment of spine metastases[END_REF], [START_REF] Deshazer | A dosimetric comparative analysis of Brainlab elements and Eclipse RapidArc for spine SBRT treatment planning[END_REF] and gradient index [START_REF] Saenz | A dosimetric analysis of a spine SBRT specific treatment planning system[END_REF] was observed. In another study comparing the same TPS, Piliero et al. [START_REF] Piliero | Plan quality and consistency in spine radiosurgery treatment planning: comparison between automatic treatment planning with Elements Spine SRS and manual inverse planning with Varian Eclipse[END_REF] concluded to a similar plan quality between both TPS. E-SRS was also compared to Monaco (Elekta LTD, Stockholm) and Pinnacle (Phillips, The Netherlands), showing a dosimetric superiority of E-SRS plans with a better ability to satisfy all the constraints simultaneously [START_REF] Saenz | A dosimetric analysis of a spine SBRT specific treatment planning system[END_REF]. Finally, the feasibility of the clinical application of the TPS was demonstrated through the evaluation of the local control in 32 patients, with an 86% local control rate at 6 and 9 months [START_REF] Giaj-Levra | Feasibility and preliminary clinical results of linac-based Stereotactic Body Radiotherapy for spinal metastases using a dedicated contouring and planning system[END_REF].

In this study, we evaluated for the first time the performances of E-SRS 3.0 TPS with comparison to plans generated with the SRS dedicated Cyberknife TPS (Multiplan 5.1.2 Accuray, Sunnyvale, USA) on a large cohort of 35 patients previously treated in our institution with the CyberKnife. Similar to other work, a comparison with the TPS that we use routinely for our standard VMAT treatments (Monaco V5.51, Elekta, Crawley, UK) was also done on the same patients' cohort [START_REF] Saenz | A dosimetric analysis of a spine SBRT specific treatment planning system[END_REF]. The objective was to investigate the real gain brought by the TPS, free from the technical benefit potentially linked to the treatment machine. The comparison of the different solutions was based on dosimetric and geometric indices. The complexity of VMAT plans provided by Monaco and E-SRS was finally assessed by calculating modulation complexity scores (MCS) [START_REF] Mcniven | A new metric for assessing IMRT modulation complexity and plan deliverability[END_REF] and analysing EPID-based quality assurance (QA) results.

Materials and methods

Patients and clinical protocol

Thirty-five patients previously treated with the CyberKnife for spine metastasis were included in this study which was approved by the hospital ethics committee. Targets and OARs were contoured according to the international guidelines [START_REF] Cox | International Spine Radiosurgery Consortium consensus guidelines for target volume definition in spinal stereotactic radiosurgery[END_REF]. The PTVs were generated by adding a 2 mm margin to the clinical target volumes (CTVs). A 2 mm margin around the OARs was added to define planning organ at risk volumes (PRV). According to UK SABR consortium [13], optimization volumes, excluding spinal cord PRV, were also generated to facilitate the optimization. In this study, the PTV coverage evaluation was performed on this PTV excluding spinal cord. The patient sample was sorted according to the vertebral level as follows: cervical (n = 10 patients), thoracic (n = 15 patients) and lumbar (n = 10 patients). The prescribed dose was 30 or 35 Gy delivered in 5 fractions to the PTV. A description of the patient cohort is provided in Table 1. The average volume of the PTVs was 38.3 ± 25.5 cm 3 . The affected spinal regions were very different from one patient to the other, leading to target volumes of varying complexity.

Treatment planning

Clinical plans were performed on MultiPlan TPS and delivered with a CyberKnife VSI robotic surgery system. Up to 3 cones from 10 mm to 50 mm diameter were used to deliver 6FFF MV multiple non-coplanar photon beams. The Monitor Units (MU) and number of nodes were optimized to limit the treatment duration from 35 to 45 min per fraction. Typically, total MU were limited to 11,000 to achieve this objective. Due to the robot geometrical limitation of this CyberKnife version, posterior beams were limited to an angle of 15° from the couch top, leading to a very limited number of posterior incidences. The final dose distribution was calculated in sequential mode using the Ray Tracing (RT) algorithm.

Each case was blindly replanned in VMAT using E-SRS 3.0 TPS. Two coplanar 360° arcs were used with collimator angle set to 80° during clockwise (CW) rotation and 280° during counterclockwise (CCW) rotation. Duplication was enabled to allow up to 4 arcs to be used to deliver the treatment, with a maximum of 90 control points (CP) per arc. The final dose distribution was calculated to medium with MC algorithm using a 2% variance and a 2 mm grid size.

VMAT plans were finally performed with Monaco 5.51 TPS using 2 CW/CCW 360° arcs, leading to 4 arcs. The collimator angles were set to 80° and 280° for each CW/CCW arc, respectively. A minimum segment width of 0.5 cm, a medium smoothing filter and a maximum of 150 CP per arc, were used, respectively, as segmentation parameters. The final dose distribution was calculated to medium with MC algorithm using a 1% variance and a 2 mm grid size.

Monaco and E-SRS treatment plans were delivered using 6 MV photon beams from a VersaHD (Elekta AB, Stockholm, Sweden) equipped with an Agility collimator with 5 mm leaves at isocenter.

Dosimetric goals

Similarly to the CyberKnife treatment plans, a dose of 30 or 35 Gy in 5 fractions was prescribed at the 80% isodose. Hence the objectives for VMAT replanning with both TPS were to respect a maximal dose that mimics a prescription on the 80% isodose. The planning objective was consistent with the UK SABR consortium guidelines for spinal metastases [13], which recommends a minimum PTV coverage with the prescription isodose of 95%. This objective was not mandatory in case of epiduritis. Hot spot of 125% of the prescribed dose was tolerated within the PTV. The dose constraints published by R. Timmerman for a 5-fractions treatment were used to set the dose constraints to the OARs, in accordance with our institution's clinical practices [START_REF] Timmerman | An overview of hypofractionation and introduction to this issue of seminars in radiation oncology[END_REF]. Hence the main constraints were to respect a maximum dose of ≤30 Gy, ≤34 Gy, ≤35 Gy and ≤38 Gy to the PRVspinal cord, PRV spinal canal, PRVesophagus and PRVtrachea, respectively, and a dose received by 200 cm 3 of kidney ≤17.5 Gy.

Dosimetric analysis

Dose distributions and Dose Volume Histograms (DVH) were analyzed in an independent software, MIM Maestro 7.0.2 (MIM software Inc., Cleveland, USA), using a dedicated workflow. All plans were resampled in MIM to a 2 mm grid size.

The different plans were compared by assessing differences in PTV V100%, D2%, D50% and D98%, and D2%, CTV D50%, and Dmax to the PRVspinal cord. Dose conformation was characterized with Paddick Conformity Index (PCI):

where TVPTV is the PTV volume covered by prescription isodose, TV is the PTV volume, and PIV is the volume of prescription isodose. Dose fall-off was evaluated using Gradient Index (GI):

where PIV50 is the volume of half the prescription isodose and PIV is the prescription isodose.

In addition, plan quality indexes (PQI) were calculated, derived from the formula proposed by Jornet et al. [START_REF] Jornet | Multicentre validation of IMRT pre-treatment verification: Comparison of in-house and external audit[END_REF], to compare planning goals achievement among TPS. Sub-PQI focusing on target metrics (PQItarget), on PRVspinal cord constraints (PQIOAR) and on the dose conformation (PQIindices) were first calculated according to the following formulas: where a value of 2 was attributed to w to give a higher importance to PTV coverage during the clinical validation of treatment plans.

For all three PQIs, goal values are reported in Table 2, and were chosen based on our clinical constraints and experience for PTV V100%, spinal cord Dmax, pCI and GI. They were derived from average values over all Cyberknife plans for each prescription scheme for PTV D98%, PTV D2% and CTV D50% assuming that Cyberknife was our gold standard in this study. Finally, a global PQI was obtained by using the following formula:

Modulation complexity score and plan quality assurance

For VMAT plans only, a modulation complexity score (MCS) was calculated from the segments' details report (segment shape, area and weight) extracted from the TPS following the formula introduced by Mc Niven et al. [START_REF] Mcniven | A new metric for assessing IMRT modulation complexity and plan deliverability[END_REF]. The number of MU and control points were also reported as indicators of modulation complexity.

To verify deliverability of the plans, QA measurements were performed for 10 patients randomly chosen, using electronic portal imaging device (EPID) and ThinkQA EPIbeam 1.0.6.21 software (DOSIsoft, Cachan, France). The passing criteria was local gamma index with 2%-2 mm tolerance with a 10% threshold. The passing rate was expected to be above 95%.

Statistical analysis

For assessing the statistical significance (p < 0.05), all differences between manual and automatically generated plans were evaluated using Wilcoxon signed ranks test.

Results

Dosimetric analysis

An example of dose distributions and dose-volume histograms obtained with the three TPS, for one patient, are presented in Fig. 1, Fig. 2, respectively, and the results for all dosimetric parameters summarized in Table 3, Table 4, Table 5 for cervical, thoracic and lumbar vertebrae, respectively. Regardless of the vertebra level, PTV V100% and PTV D98% were not statistically different between all the TPS (p > 0.107), except for the thoracic level where Monaco showed better PTV V100% than Multiplan (p = 0.026) (Table 4). Conversely, PTV D50% of E-SRS treatment plans were found to be significantly higher (p < 0.004) in 35 and in 34 out of 35 treatment plans, compared to Multiplan and Monaco, respectively, with PTV D50% being up to 1.8 Gy higher than with other solutions at the lumbar level. On the contrary, PTV D50% were not significantly different between Monaco and Multiplan (p > 0.063), with the exception of the thoracic level, where it was significantly superior with Monaco than with Multiplan (p = 0.015). CTV D50% were also significantly higher with E-SRS TPS compared to other TPS (p < 0.031), whatever the vertebra level, with D50% to the CTV being up to 1.6 Gy higher than with other TPS. In contrast, CTV D50% were not significantly different between Monaco and Multiplan (p-value > 0.063). Consequently, PTV D2% was found significantly higher with E-SRS than with other TPS (p < 0.004) and also significantly different between Monaco and Multiplan (p < 0.002) with a PTV D2% slightly lower in Monaco than in Multiplan, with the exception of the cervical level (p = 0.131). Regarding the PCI, no significant difference was observed between the 3 vertebrae levels. On average, at equivalent PTV coverage, higher PCI indices were obtained at the lumbar level for all the TPS. The PCI were not statistically different for all vertebral levels, with the exception of the lumbar level, where E-SRS plans resulted in better PCI compared to Monaco (p = 0.032) and Multiplan (p = 0.024).

E-SRS plans also resulted in better GI than Monaco plans, whatever the vertebral level (p < 0.009), and better GI than Multiplan plans for the thoracic level only (p = 0.007) (Table 4). Finally, the GI was significantly better in Multiplan than in Monaco for the cervical level only (p = 0.006), whereas no statistical difference was observed on other vertebral levels (p ≥ 0.232).

The Dmax to the spinal cord was respecting the dose constraints criteria for all the treatment plans (Table 3, Table 4, Table 5), and no significant difference was found between the TPS, except between Multiplan and E-SRS on the thoracic level (p = 0.008). In contrast, the D2% to the spinal cord was found to be significantly lower with E-SRS compared to Multiplan in 30 out of 35 treatment plans (p < 0.037), and compared to Monaco at the thoracic level only, for 11/15 instances (p = 0.038). Regarding Monaco and Multiplan TPS, D2% to the spinal cord was significantly lower for Monaco compared to Multiplan for the cervical level (p = 0.037), while there was no significant difference observed for other vertebrae levels (p > 0.188).

Regarding the number of MU, on average, Multiplan plans had up to 4.5 times and 3 times more MU than E-SRS and Monaco plans, respectively (Table 6). Finally, the PQIglobal scoring did not reflect a significant superiority of Cyberknife plans over E-SRS and Monaco plans, with a better scoring obtained on 13 patients versus 12 and 10 patients for E-SRS and Monaco, respectively (Fig. 3).

Furthermore, for several patients (i.e. P5, P10, P12, P13 and P33), the PQIglobal of Cyberknife plans were significantly lower than for VMAT plans. The same observation was done for some Monaco plans (i.e. P8, P11, P27) in comparison with the two other techniques. Note that for P24, Cyberknife and Monaco plans had a significantly lower score than E-SRS plans. Finally, E-SRS plans presented the narrowest dispersion of points over all the patients' cohort.

PQItarget scoring reflected the superiority of Monaco plans compared with E-SRS and Cyberknife, with 22 patients over 35 having better scoring than the two other modalities. In contrast, the PQIindices calculation showed the superiority of the 2 other techniques compared with Monaco, each of the 2 techniques having obtained the best score for 11 patients. Finally, PQIoar scoring showed the superiority of E-SRS and Monaco plans over Cyberknife plans, with a better scoring obtained on 16 and 15 patients for E-SRS and Monaco, respectively, versus 4 for Cyberknife plans (Fig. 3).

QA results

VMAT E-SRS and Monaco plans QA results are shown in Table 7. Despite a number of MU up to 1.7 times higher with Monaco than with Elements Spine SRS, the MCS were found not to be statistically different. However, E-SRS plans generally resulted in lower MCS. Consequently, even if not statistically significant, the mean gamma passing rate was slightly lower for E-SRS (95.8%) compared with Monaco TPS (96.7%) using a 2%-2mm local analysis method. For Monaco TPS one single plan showed a gamma index < 95% (91.3%). Similarly, for Elements Spine SRS, for 3 patients, the gamma passing rate was comprised between 91.2 and 94.3%. However, all QA plans would have been validated in a clinical situation with global 3%-3mm gamma.

Discussion

The main purpose of this study was to investigate for the first time the performances of a dedicated TPS performing SRS treatment plans delivered with a conventional Linac, compared to a dedicated Linac for the management of SRS, the Cyberknife and its TPS. The comparison was based on a cohort of 35 patients with spine tumours previously treated with the Cyberknife, in our institution. A comparison with Monaco TPS that we use routinely for our standard VMAT treatment was also performed.

Using VMAT delivery technique, E-SRS and Monaco TPS produced high quality plans with an achievement of our dosimetric goals comparable to the CyberKnife Multiplan TPS that was our "gold standard" modality. All plans were matching the dosimetric goals. The PTV coverage, Dmax to the spinal cord, and conformation indices were generally found to be comparable between all modalities and the vertebral level did not have a significant impact on dosimetric results.

There was no significant difference between Monaco and Multiplan TPS treatment plans, with the exception of the PTV D50% to the thoracic level, which was higher with Monaco. Hence, with Monaco TPS, we were able to mimic the quality of the plans obtained with Multiplan TPS. This was made possible because the prescription templates used in Monaco were based on the analysis of the dosimetric criteria collected from a large cohort of patients (>100) previously treated with the CyberKnife. Briefly, our prescription template was developed by considering a dose prescribed at the 80% isodose. The maximum dose allowed in the target was both considered as a constraint and as an objective. Finally, the prescribed dose to the CTV was generally 1.11 to 1.13 times higher than the dose prescribed at the 80% isodose to match the D50% generally obtained with the CyberKnife. Thus, although Monaco TPS is a VMAT dedicated TPS, initially designed to plan with homogeneity requirements specific to this technique, we showed that it was possible to develop prescriptions specific to SBRT which give similar dosimetric results to the CyberKnife technique [START_REF]ICRU Report 83, Prescribing, Recording, and Reporting Intensity-Modulated Photon-Beam Therapy (IMRT) -ICRU n[END_REF], [START_REF]ICRU Report 91, Prescribing, Recording, and Reporting of Stereotactic Treatments with Small Photon Beams -ICRU n[END_REF].

Conversely, the D50% to the PTVs and CTVs were found significantly higher (up to 1.8 Gy higher) for E-SRS compared to Multiplan and Monaco treatment plans. The D2% to the spinal cord was significantly lower with E-SRS than with Multiplan and Monaco, leading to a spinal cord DVH clearly favourable regarding the dose-volume constraints proposed in the literature [13], [START_REF] Timmerman | An overview of hypofractionation and introduction to this issue of seminars in radiation oncology[END_REF]. Thus, in the optimization algorithm, a harder constraint was probably applied to the maximum dose to the spinal cord, compared with other TPS. This stricter competition between PTV and OAR combined with the use of conformation and gradient objectives in the prescription, may have played an important role in the increase of median dose to PTV and CTV observed on E-SRS plans. Note that the presence of these parameters did not result in statistically significant differences between TPS on the calculated PCI and GI indices.

The PQI calculations results were consistent with above observations. Although several dosimetric goals came from averaged Cyberknife's values obtained over all patients' cohorts, the PQIglobal scoring did not highlight a superiority of Cyberknife compared with both VMAT plans. Note that the goals and hypothesis made for PQI calculations may strongly impact the trends observed and are very subjective and dependant on the relative importance given to a constraint or to another, related to clinical practices. In our study, a weighting factor was affected to PTV V100% only, as this criterion has the most relative importance in the clinical decision in our practice. In addition, PQIindices had the same relative importance than PQIoar and PQItarget in our PQIglobal calculations, as we considered pCi and GI metrics of importance for stereotactic plan evaluation. In an interesting way, the dispersion of the different PQIs observed for E-SRS plans was systematically smaller, suggesting a more reproducible planning strategy provided by automatization.

Several studies have compared the CyberKnife with VMAT for Spine SBRT planning, but, to our knowledge, this is the first time that a comparison is made with E-SRS TPS [START_REF] Nalichowski | Single fraction radiosurgery/stereotactic body radiation therapy (SBRT) for spine metastasis: A dosimetric comparison of multiple delivery platforms[END_REF], [START_REF] Moustakis | Treatment planning for spinal radiosurgery: A competitive multiplatform benchmark challenge[END_REF], [START_REF] Esposito | SBRT planning for spinal metastasis: indications from a large multicentric study Strahlenther[END_REF]. One study was performed on an anthropomorphic phantom and compared the Cyberknife technology to the Tomotherapy (Accuray, Sunnyvale, USA), Vero (Brainlab, Munich, Germany) and VMAT planning using Eclipse (Varian, Palo Alto, USA) and concluded that the Cyberknife had the sharpest dose falloff and achieved better plan quality for small targets [START_REF] Nalichowski | Single fraction radiosurgery/stereotactic body radiation therapy (SBRT) for spine metastasis: A dosimetric comparison of multiple delivery platforms[END_REF]. In a multicentric study, 3 different patients with various spinal lesions were planned using the Tomotherapy and different VMAT planning solutions and compared to previously published Cyberknife benchmark data [START_REF] Moustakis | Treatment planning for spinal radiosurgery: A competitive multiplatform benchmark challenge[END_REF]. They concluded that VMAT plans were able to outperform these established Cyberknife benchmark data [START_REF] Blanck | Inverse treatment planning for spinal robotic radiosurgery: an international multi-institutional benchmark trial[END_REF].

Several studies have compared E-SRS with other VMAT solutions [START_REF] Trager | Evaluation of elements spine SRS plan quality for SRS and SBRT treatment of spine metastases[END_REF], [START_REF] Deshazer | A dosimetric comparative analysis of Brainlab elements and Eclipse RapidArc for spine SBRT treatment planning[END_REF], [START_REF] Piliero | Plan quality and consistency in spine radiosurgery treatment planning: comparison between automatic treatment planning with Elements Spine SRS and manual inverse planning with Varian Eclipse[END_REF], [START_REF] Saenz | A dosimetric analysis of a spine SBRT specific treatment planning system[END_REF]. In a first study, similar dosimetric results were found on the dose to the spinal cord and to the PTV, when comparing the E-SRS with Eclipse TPS (AAA algorithm) [START_REF] Trager | Evaluation of elements spine SRS plan quality for SRS and SBRT treatment of spine metastases[END_REF]. The authors also noted an improvement of PCI and GI when using Elements Spine SRS, with the same beam geometry in both TPS but with duplication option in Spine SRS. In another study, the same TPS were compared, but using different algorithms for final dose calculation (AcurosXB for Eclipse TPS; Pencil Beam and MC algorithms for Spine SRS) [START_REF] Deshazer | A dosimetric comparative analysis of Brainlab elements and Eclipse RapidArc for spine SBRT treatment planning[END_REF]. Similar to us, they found a significant difference on spinal cord sparing, and on the maximum dose to the CTV with on average a maximum dose 1.6 Gy higher in Spine SRS using MC algorithm. Again, the same beam geometry was used in both TPS. However, contrary to us, similar PCI was found between their TPS used in routine for VMAT and both Element's algorithms. In another study, authors also compared Eclipse plans calculated with AcurosXB with E-SRS plans calculated with MC algorithm [START_REF] Piliero | Plan quality and consistency in spine radiosurgery treatment planning: comparison between automatic treatment planning with Elements Spine SRS and manual inverse planning with Varian Eclipse[END_REF]. One major methodological difference between this study and others was that planners were free to choose their ballistic in Eclipse (2 complete arcs with different collimator geometry), while E-SRS plans were computed with a defined geometry limited to one arc with duplication on. Similar to us and previous author's results, they found median doses to target volume 0.8 to 1.5 Gy higher in Spine SRS than in Eclipse, but they did not find a better spinal cord sparing with E-SRS and they conclude to an overall comparable plan quality between both TPS. Finally, Saenz et al. [START_REF] Saenz | A dosimetric analysis of a spine SBRT specific treatment planning system[END_REF] compared E-SRS (MC) plans to Pinnacle (adaptive convolve algorithm) and Monaco plans (MC) and they found all their dosimetric criteria to be better achieved with Elements Spine SRS. In this study, target high dose was characterized by PTV D5% and the three TPS were respecting the constraint of D5% < 25 Gy for 20 Gy prescribed dose. Contrary to other studies, the increase in the median dose to the targets observed with E-SRS was not observed here.

Thus, the increase in the median dose to PTV and CTV and spinal cord sparing in Spine SRS was a shared observation between most authors. However, we do not have any objective judgment about whether this dose increase to CTV and PTV should be clinically desirable or not, in a stereotactic prescription scheme, where maximal and median doses to the target volumes can vary considerably between centers using different prescription strategies. Spinal cord sparing is also a debatable question as the gradient augmentation has to be related to the patient repositioning and monitoring strategy and may expose to a higher risk of recurrence, even if an acute spinal cord sparing allows easier re-irradiation. In our institution, spinal cord sparing is not pushed as far as in Spine SRS, even with CyberKnife.

One of the factors that can influence the results is the inter-operator variability. As shown in a multi-institutional study focused on CyberKnife planning, a wide variability between treatment planning approaches and dosimetric results was observed between participants [START_REF] Blanck | Inverse treatment planning for spinal robotic radiosurgery: an international multi-institutional benchmark trial[END_REF]. The plan quality dependence on planners' skill was also discussed elsewhere [START_REF] Esposito | SBRT planning for spinal metastasis: indications from a large multicentric study Strahlenther[END_REF] and may have more impact than the TPS and treatment device combination. In this latter study, the authors proposed a crowd knowledge-based approach in spine stereotactic planning and showed this method to be effective in improving the global plan quality. To reduce plan quality dependence to user's skill, automated planning may be another promising solution.

As observed in our study, the PQI results were more homogeneous with the dedicated TPS than with manual planning solutions. The E-SRS is making a step in more reproducible plan solutions, integrating semi-automated planning, with a very limited number of parameters available for the operator, during the optimization process. The advantage of this TPS developed exclusively for SRS planning lies in its user-friendliness and its optimization of all phases, including: time-consuming delineation tasks simplified by anatomical libraries integrated into the software, allowing this phase to be automated [START_REF] Rogé | Evaluation of a dedicated software "Elements™ Spine SRS, Brainlab®" for target volume definition in the treatment of spinal bone metastases with stereotactic body radiotherapy[END_REF], [START_REF] Chen | Clinical evaluation of an auto-segmentation tool for spine SBRT treatment[END_REF]; optimization constraints automatically managed once they are set by the user for each fractionation protocol. Furthermore, the optimization duplication option, which lead to the target subdivision into simpler volumes to avoid concavities/convexities, seems to allow a high level of spinal cord preservation.

When performing dosimetric studies, a particular attention should be paid to the tools used for the dose distribution comparison. In this study, a third-part software, MIM, with a systematic resampling in a 2 mm grid size, was used as unique tool to compare dose distributions. This methodology can lead to slight differences in DVH values observed in TPS and in MIM. In our results, this was particularly true for E-SRS plans, which have an adaptive sub-millimetric calculation grid, with a resampling of the grid depending on the size of the structure. Users should be aware of finer and non-isomorphic grids of Elements TPS. Using a unique software for dose distribution comparison allows to compare plans with the same tools, methods and resolution. Indeed, a study of Stanley et al. found differences in reported metrics for Elements, Eclipse, MIM, Raystation and Velocity when importing Dicom RTDose and RTStruct files in these different software [START_REF] Stanley | Accuracy of dose-volume metric calculation for small-volume radiosurgery targets[END_REF]. Significant variations in volume calculations were also observed elsewhere [START_REF] Ma | Reliability of contour-based volume calculation for radiosurgery[END_REF]. Resampling errors, which potentially occur when contour points from one treatment planning system are imported into another, can significantly affect volume calculations for small objects, and thus lead to discrepancies in DVH calculation and plan comparison.

Another element to consider in this study was that the comparison between dose distributions was computed with different algorithm classes. Indeed, Monaco and E-SRS dose distribution were computed with a MC algorithm, whereas Multiplan dose distribution was computed with the RT algorithm. It has been demonstrated that MC algorithms provide a higher accuracy in dose calculation in heterogeneous medium and are recommended for stereotactic treatment planning [START_REF]ICRU Report 91, Prescribing, Recording, and Reporting of Stereotactic Treatments with Small Photon Beams -ICRU n[END_REF]. Hence, in two previous studies, RT and MC algorithms were compared for Cyberknife' spine planning. Li et al. observed differences in PTV coverage from 1% in cervical region up to 13% in thoracic region [START_REF] Li | Yang Comparison of ray tracing and Monte Carlo calculation algorithms for spine lesions treated with cyberknife[END_REF]. The dose to the spinal cord found with MC algorithm was generally lower than with RT algorithm. Similarly, Okoye et al. found a significant loss of PTV coverage when recalculating in MC thoracic spinal plans optimized in RT [START_REF] Okoye | Comparison of ray tracing and Monte Carlo calculation algorithms for thoracic spine lesions treated with CyberKnife-based stereotactic body radiation therapy[END_REF]. Finally, both authors concluded in the necessity to use MC algorithm for thoracic spinal lesions planning. In our study, MC was not clinically implemented at the time the plans were performed. The goal of this study was not to evaluate the accuracy of the algorithms but to compare the optimization results produced by three TPS, and to assess the E-SRS automatic planning performances. However, a study conducted in parallel allowed us to obtain the same trends mentioned by the two previous studies. In this study, the Cyberknife plans were retrospectively recomputed with an independent MC dose engine (Surecalc v. 7.3.1; scimoca 1.7.1.5050 dose engine (MIM software Inc., Cleveland, USA)). We observed that PTV coverage and OAR Dmax were overestimated with the RT algorithm, whatever the vertebral level (Fig. 1 -Supplementary Material). PTV V100% was found to be on average 10%, 15% and 5% lower for cervical, thoracic and lumbar vertebral level, respectively, with MC calculations compared with RT. D2% to spinal cord followed the same trend, with an average dose reduction of up to 8% to the thoracic level.

In the present study, we observed that E-SRS plans resulted in a significantly lower number of MU and in a lower MCS compared to Monaco TPS but without compromising the plan delivery quality. Compared with Multiplan TPS, the number of MU obtained for the plans performed with E-SRS and Monaco TPS were significantly lower, which resulted in an improved delivery efficiency for an equivalent plan quality compared with our gold standard. This also resulted in an important decrease of the beam on time for the VMAT technique compared with the Cyberknife. This was estimated to 6 min, 10 min and 40 min, for Elements Spine SRS, Monaco and Cyberknife treatments, respectively.

Finally, the CyberKnife clinical plans were performed using fixed cones. Since then, our CyberKnife was changed for a S7 model implying changes in our planning method and delivery for spine SBRT, i.e the use of VOLO optimization algorithm (Precision, Accuray, Sunnyvale) and the IRIS as beam collimation. These recent features added to Cyberknife linac and TPS optimize the beam arrangement, decrease the beam on time and the number of MU [START_REF] Zeverino | Novel inverse planning optimization algorithm for robotic radiosurgery: First clinical implementation and dosimetric evaluation[END_REF], [START_REF] Schüler | Clinical impact of the VOLO optimizer on treatment plan quality and clinical treatment efficiency for CyberKnife[END_REF]. Furthemore the use of the MLCincise could contribute to a treatment delivery time decrease of 30% compared to fixed cones [START_REF] Kim | Clinical outcomes of multileaf collimator-based CyberKnife for spine stereotactic body radiation therapy[END_REF]. However, our goal in this study was to assess E-SRS performances for spinal SBRT automatic planning compared to plans produced by skilled users on other TPS, based on our clinical background.

Conclusion

Spinal SBRT planning with E-SRS offers very qualitative dosimetric solutions, close to our clinical experience with both Multiplan and Monaco TPS. There is a slightly different balance between goals to spinal cord sparing and PTV doses between our traditional way of prescribing in our clinical TPS and the natural Spine SRS optimizer behaviour but, to date, there is no any objective judgment about whether these differences are clinically desirable or not, in a stereotactic prescription scheme. Finally, the E-SRS special approach offers a very effective and user-friendly semi-automated planning tools and is secure and promising for linac-based spinal SBRT.

Fig. 1 .Fig. 2 .

 12 Fig. 1. Axial, sagittal and coronal views of Monaco, Multiplan and Elements Spine SRS dose distributions for one patient. The CTV is in red, the PTV is in blue, and the spinal canal in cyan. 100% isodose line is thick orange line, 50% isodose line is thick cyan line (prescribed dose: 30 Gy in 5 fractions) and 125% isodose line is in thick red. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 3 .

 3 Fig. 3. PQI calculated for global plan quality scoring (A), target constraints achievement scoring (B), indices scoring (C) and OAR constraints achievement scoring.

Table 1 .

 1 Cohort description.

	Number of patients

PTV volume (cc) mean ± standard deviation (interval)

	Cervical	16.17 ± 7.39 (6.64-32.35)
	Thoracic	39.52 ± 25.06 (13.08-99.43)
	Lumbar	58.65 ± 20.6 (43.06-100.98)
	Dose de prescription
	30 Gy in 5 fractions	23 (65.71%)
	35 Gy in 5 fractions	12 (34.29%)

Table 2 .

 2 Goals list and weighting factor for the metrics involved in PQI calculations.

	Structures/Indices	Goal

Table 3 .

 3 Summary of dosimetric results for cervical vertebrae plans. Values are average across the ten cases.

	Cervical Vertebrae Plans		E-SRS	Multiplan	Monaco
			Mean [range]	Mean [range]	Mean [range]
	PTV V100% (%)		96.1 [93.0-98.5]	94.5 [85.3-99.6]	96.2 [92.5-98.9]
	30 Gy (8 patients)	PTV D2% (Gy)	37.3 [37.2-37.4]	36.6 [36.1-37.4]	36.2 [35.4-37.3]
		PTV D50% (Gy)	35.0 [34.3-35.6]	33.5 [32.4-34.3]	33.8 [33.1-34.3]
		PTV D98% (Gy)	29.4 [28.6-30.3]	28.6 [20.4-30.6]	28.6 [20.1-30.4]
		CTV D50% (Gy)	35.7 [35.1-36.2]	34.1 [32.9-35.6]	34.4 [33.7-35.0]
	35 Gy (2 patients)	PTV D2% (Gy)	43.3 [43.3-43.4]	42.9 [42.9-43.0]	42.7 [42.6-42.7]
		PTV D50% (Gy)	39.9 [39.8-40.2]	39.1 [38.6-39.5]	39.5 [38.9-40.0]
		PTV D98% (Gy)	33.6 [33.3-33.9]	32.5 [31.5-33.4]	33.9 [33.6-34.1]
		CTV D50% (Gy)	41.0 [40.9-41.2]	40.2 [40.1-40.2]	40.3 [39. 8-40.8]
		CI	0.79 [0.73-0.87]	0.77 [0.68-0.86]	0.77 [0.58-0.84]
		GI	3.95 [3.60-4.29]	3.97 [3.13-5.08]	5.13 [4.07-6.65]
	Dose to Spinal Cord PRV	Dmax (Gy)	29.0 [26.6-31.3]	29.3 [25.9-32.6]	28.9 [25.6-31.8]
		D2% (Gy)	23.9 [20.2-26.9]	26.2 [23.3-28.6]	24.5 [21.3-28.8]

Table 4 .

 4 Summary of dosimetric results for thoracic vertebrae plans. Values are average across the fifteen cases.

	Thoracic Vertebrae Plans		E-SRS	Multiplan	Monaco
			Mean [range]	Mean [range]	Mean [range]
	PTV V100% (%)		96.1 [92.0-98.8]	95.1 [89.6-98.4]	96.7 [91.0-98.7]
	30 Gy (11 Patients)	PTV D2% (Gy)	37.2 [36.3-37.5]	36.6 [36.4-36.8]	36.3 [35.3-36.7]
		PTV D50% (Gy)	35.2 [34.5-35.6]	33.5 [32.8-34.3]	34.0 [33.2-34.5]
		PTV D98% (Gy)	29.4 [28.7-30.6]	29.5 [28.2-30.1]	29.7 [28.7-30.3]
		CTV D50% (Gy)	35.6 [35.0-35.9]	34.2 [33.5-34.9]	34.4 [33.5-35.1]
	35 Gy (4 patients)	PTV D2% (Gy)	43.4 [43.3-43.4]	42.8 [42.6-43.1]	42.1 [40.9-42.7]
		PTV D50% (Gy)	40.9 [39.8-41.5]	39.3 [39.1-39.5]	39.5 [38.6-40.2]
		PTV D98% (Gy)	33.1 [32.4-34.6]	33.0 [32.5-33.3]	33.3 [30.2-34.7]
		CTV D50% (Gy)	41.4 [40.6-41.9]	40.1 [39.8-40.4]	40.0 [38.8-40.8]
		CI	0.80 [0.67-0.89]	0.78 [0.69-0.84]	0.78 [0.61-0.88]
		GI	3.75 [2.99-4.43]	4.31 [3.47-5.66]	4.22 [3.45-5.48]
	Dose to Spinal Cord PRV	Dmax (Gy)	29.4 [28.7-30.8]	29.9 [28.6-30.6]	29.4 [26.1-30.8]
		D2% (Gy)	25.4 [19.7-27.3]	27.1 [20.0-27.7]	26.4 [21.2-29.8]

Table 5 .

 5 Summary of dosimetric results for lumbar vertebrae plans. Values are average across the ten cases.

	Lumbar Vertebrae Plans		E-SRS	Multiplan	Monaco
			Mean [range]	Mean [range]	Mean [range]
	PTV V100% (%)		96.0 [94.7-98.4]	96.0 [90.6-99.8]	96.6 [94.8-99.0]
	30 Gy	PTV D2% (Gy)	37.34 [37.1-37.5]	36.59 [36.3-36.8]	36.2 [35.6-36.5]
	(4 patients)
		PTV D50% (Gy)	35.39 [34.8-35.8]	33.91 [33.1-35.2]	33.9 [33.4-34.2]
		PTV D98% (Gy)	29.35 [28.5-30.3]	30.20 [29.6-30.7]	29.8 [28.4-30.4]
		CTV D50% (Gy)	35.59 [35.3-35.9]	34.01 [33.6-34.4]	34.5 [34.5-34.6]
	35 Gy	PTV D2% (Gy)	43.54 [43.3-43.8]	42.68 [42.4-42.9]	42.3 [42.0-42.6]
	(6 patients)
		PTV D50% (Gy)	41.16 [40.5-41.5]	39.38 [38.7-39.8]	39.9 [39.7-40.2]
		PTV D98% (Gy)	33.21 [31.7-34.1]	33.30 [31.8-34.6]	32.8 [31.4-33.4]
		CTV D50% (Gy)	41.61 [40.7-42.1]	39.8 [39.3-40.5]	40.1 [39.3-40.7]
		CI	0.83 [0.77-0.89]	0.79 [0.71-0.84]	0.83 [0.76-0.92]
		GI	3.60 [3.38-3.96]	3.75 [3.26-4.90]	3.97 [3.60-4.51]
	Dose to Spinal Cord PRV	Dmax (Gy)	29.1 [27.9-31.6]	29.8 [27.8-31.8]	28.9 [26.3-30.9]
		D2% (Gy)	25.0 [22.8-26.6]	26.6 [22.6-29.3]	25.7 [22.9-28.4]

Table 6 .

 6 Summary of average total MU for the three TPS, per vertebral level and for all patients.

	Empty Cell	E-SRS	Multiplan	Monaco
		Total MU: Mean [range]	Total MU: Mean [range]	Total MU: Mean [range]
	Cervical plans	7476 [5730-12385]	30,850 [18242-43251]	10,365 [7189-15050]
	Thoracic plans	9233 [6950-12880]	41,546 [16303-76623]	14,614 [9939-22092]
	Lumbar plans	9845 [7620-18700]	46,998 [35376-56432]	16,495 [10246-21705]
	All plans	8906	40,048	13,938

Table 7 .

 7 Summary of QA and modulation complexity scores (MCS) results for 10 patients.

	Analysis Criteria	E-SRS	Monaco	p-value
		Mean [Range]	Mean [Range]
	Mean Local Gamma Index	0.46 [0.41-0.53]	0.42 [0.35-0.49]	0.074
	2%-2 mm
	Mean Global Gamma Index	0.39 [0.35-0.45]	0.35 [0.29-0.43]	0.123
	2% -2 mm
	Local 2% -2 mm Gamma Index passing rate (%)	95.8 [91.2-98.3]	96.7 [91.3-99.3]	0.232
	Global 2%-2mm Gamma Index passing rate (%)	97.4 [92.7-99.4]	97.9 [92.6-99.4]	0.770
	MCS	0.25 [0.14-0.37]	0.27 [0.20-0.34]	0.432
	MU	8323 [5775-12385]	12,008 [8017-20184]	0.002

Acknowledgements

This work was performed within the framework of the SIRIC LYriCAN Grant INCa-INSERM-DGOS-12563, and the LABEX PRIMES (ANR11-LABX-0063) of Université de Lyon, within the program Investissements d'Avenir (ANR-11-IDEX-0007) operated by the ANR. We are grateful to Rebecca Ljungqvist and Mathilde Favier from Brainlab society and to Stefano Sobkowiak from MIM software society, and Sophie King for their support and advices.

Declaration of Competing Interest

The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: This work was performed in the framework of a research cooperation agreement with Brainlab society.

Appendix A. Supplementary data


Samy Gauthier

Alexy Tran-Dinh

Morilla
email: ian.morilla@uma.es

Methods Plasma proteome dynamics of COVID-19 severity learnt by a graph convolutional network of multi-scale topology

Efforts to understand the molecular mechanisms of COVID-19 have led to the identification of ACE2 as the main receptor for the SARS-CoV-2 spike protein on cell surfaces. However, there are still important questions about the role of other proteins in disease progression. To address these questions, we modelled the plasma proteome of 384 COVID-19 patients using protein level measurements taken at three different times and incorporating comprehensive clinical evaluation data collected 28 d after hospitalisation. Our analysis can accurately assess the severity of the illness using a metric based on WHO scores. By using topological vectorisation, we identified proteins that vary most in expression based on disease severity, and then utilised these findings to construct a graph convolutional network. This dynamic model allows us to learn the molecular interactions between these proteins, providing a tool to determine the severity of a COVID-19 infection at an early stage and identify potential pharmacological treatments by studying the dynamic interactions between the most relevant proteins.

Introduction

The sudden spread of severe acute respiratory syndrome coronavirus 2 (Sars-Cov-2) worldwide meant one of the most remarkable public health crises of recent times. Currently, more than 60 million individuals have been infected, causing over 1.5 million deaths related with severe complications of the COVID-19 disease.

There exists proved evidence on how certain proteins such as ACE2 receptor or TMPRSS2 are used by Sars-Cov-2 as entrance gates to infect the cell via membrane fusion and endocytosis [START_REF] Ou | Characterization of spike glycoprotein of SARS-CoV-2 on virus entry and its immune cross-reactivity with SARS-CoV[END_REF][START_REF] Zhou | A pneumonia outbreak associated with a new coronavirus of probable bat origin[END_REF]. Likewise, there are also multiple clues on a likely participation of other proteins in the downstream of the disease during its progression [START_REF] Blanco | In silico mutagenesis of human ACE2 with S protein and translational efficiency explain SARS-CoV-2 infectivity in different species[END_REF][START_REF] Yang | Molecular interaction and inhibition of SARS-CoV-2 binding to the ACE2 receptor[END_REF][START_REF] Cuervo | ACE2: Evidence of role as entry receptor for SARS-CoV-2 and implications in comorbidities[END_REF][START_REF] Scudellari | How the coronavirus infects cells -and why Delta is so dangerous[END_REF]. All these experimental efforts aim to characterise the progression of COVID-19 from an in situ baseline analysis of proteomic profiles. An approach that is getting more popular nowadays is to switch off those protein interactions resulting essential to the viral infection. Thus, the targeting of protein-protein interaction interfaces may be used to discover anti-COVID-19 treatment [START_REF] Xiu | Inhibitors of SARS-CoV-2 entry: Current and future opportunities[END_REF][START_REF] Yang | Molecular interaction and inhibition of SARS-CoV-2 binding to the ACE2 receptor[END_REF]. Unfortunately, this and most of those analyses tend to overlooking nonlinear programmes of interaction determined by subsets of proteins already described in such studies. Some of those programmes are merely contributing to an innocuous reconfiguration of the secondary immune response system, but others can be causally provoking a worsening in the severity of the symptoms (i.e., hyperinflammatory syndrome) during the disease progression.

In this work, we learn the latter through graph convolutional networks (GCNs) calibrated with higher topological features extracted from raw plasma proteomic data collected by the Massachusetts General Hospital Emergency Department COVID-19 (https://www.olink.com/mgh-covid-study/- [START_REF] Filbin | Data provided by the MGH Emergency Department COVID-19 with Olink Proteomic[END_REF]). The examination of plasma as a potential source of insight into the evasion mechanisms of severe acute respiratory SARS-CoV-2 has been the subject of recent research [START_REF] Cabrera-Garcia | Plasma biomarkers associated with survival and thrombosis in hospitalized COVID-19 patients[END_REF][START_REF] Zhao | Plasma biomarkers for systemic inflammation in COVID-19 survivors[END_REF]. Plasma, as the liquid component of blood, is a complex mixture of substances, including antibodies and proteins produced by the host immune system in response to a stimulus, such as infection or disease. Alterations in the levels of these substances in the plasma may provide evidence of underlying conditions or disorders, even in the absence of direct evaluation of immune system cells.

Thus, we clustered and identified patient phenotypes by the World Health Organization (WHO)-mediated scores along their comorbidity effectors (if any) in a non-supervised fashion first and computed later with persistent homology of more than 1; 400 protein profiles in blood enhanced in endothelial cells (see Fig 1) across samples [START_REF] Xu | Finding cosmic voids and filament loops using topological data analysis[END_REF]. To the best of our knowledge, this study is the first to chart longitudinal associations between plasma protein interactions and disease outcomes in patients with COVID-19 disease. Finally, we are convinced our results will be instrumental in a later experimental validation.

Results

Study design

Our models achieved tracking protein interactions occurring postinfection by which different levels of severity developed by 384 individuals suffering from COVID-19 symptoms might be explained [START_REF] Filbin | Data provided by the MGH Emergency Department COVID-19 with Olink Proteomic[END_REF]. In this sense, the COVID status of inpatients was tested positive prior to enrolment or during hospitalisation. Then, based on that test, we discriminated 306 patients as Covid + and 78 patients as Covid -(see Fig 2A). In addition, we wanted to improve our models' interpretability by means of the detailed clinical outcomes available from each patient at day 28 of their stay in the hospital. In total, those variables encompass up to 40 different types of multi-variated sequences (see Supplemental Data 1). Those data allowed us to compute a precise overall score of severity based on discrete WHO scores (Organisation world health, 2021) provided in the cohort for each patient over time of stay in the hospital. First, we vectorised those scores with the aim of being using entropy [START_REF] Gray | Entropy and Information Theory[END_REF] to calculate the information encapsulated by the WHO scores in each patient. That information in bits [START_REF] Murphy | Machine Learning: A Probabilistic Perspective (Adaptive Computation and Machine Learning Series[END_REF] enabled the construction of a probability density function that basically well stratified individuals according to severity progression of the disease which ultimately reinforced the explanatory power of our learning omic models. Hierarchical models of WHO scale-based entropy information stratify patients by severity WHO provides scores that monitor disease progression during a patient's stay at a hospital. These scores, recorded using the discrete measure W s = f1; 2; …; 6g 2 Z + at days 0, 3, 7, and 28 (see Table 2 and Supplemental Data 1), provide a local snapshot of an individual's disease progression (see Fig 2B). To better understand and track disease progression, we propose using information theory scores, as proposed by Shannon in his seminal paper [START_REF] Shannon | A mathematical theory of communication[END_REF]. This approach allows us to efficiently store and retrieve relevant information and create a continuous distribution to model disease progression. By vectorising the WHO scores and calculating entropy, we can extract self-information and fit the transformed data (see Fig S1A ) to the best distribution by checking a comprehensive set of probability density functions. In particular, we obtain a discrete random variable V, with possible outcomes v 1 ; …; v 6 , which occur with probability Pðv 1 Þ; :::; Pðv 6 Þ and computed the entropy of each individual per vector by applying the formula 2D) and was used to calculate the bits of information [START_REF] Mackay | Information Theory, Inference, and Learning Algorithms[END_REF] for a context channel (i.e., pðv v 9 Þ) with v and v 9 in the discrete alphabet W s . Thus, we can accurately track disease progression and identify potential trends.

Hð V Þ = -å 6 i = 1 P vi Þlog bits P vi Þ À À (see Fig 2C).
The generalised continuous normal random variable of subdomain D = ½1; 2 2 R + was calculated as follows: S r ðA;

IÞ i = å N l = 1 1 D 1=2 k;k A k;l 1 D 1=2 l;l I l ; fðx; βÞ = β 2Γ 1 β exp -kxk β ;
where β 2 ½0:01; 0:99 and Γ is the function gamma [START_REF] Sun | A note on the generalized gamma function[END_REF]. By means of this function, we calibrated hierarchical clustering models that were computed by applying the Hdbscan algorithm [START_REF] Campello | Density-based clustering based on hierarchical density estimates[END_REF] From this stratification, 13 out of 384 patients were excluded to be considered as outliers with noisy data. These patients displayed dissimilar symptoms and unmatched characteristics amongst them to be included in any of the groups (see Video 4 and Fig S2).

The latent space of clinical features explains patients' stratification

We unified many local perspectives of the clinical dataset to explain models of severity These plots display inflection points in predicted age and hypertension as $Age cat$ and HTN (oldness by years on average and hypertension complaint per individual in the cohort) changes. Vertical dispersion at a single category of Age (resp. HTN) represents interaction effects with other features. To help reveal these interactions, we coloured by Fever (resp. BMI). We passed the whole explanation tensor to the colour argument in the dependence plots to pick the best feature to colour by. In this case, it selected fever symptoms (resp. Body Mass Index) because that highlights that the average age (hypertension) per inpatient has more (less) impact on Covid severity for categories with a low (high) Fever (BMI) value.

the predicted risk of experiencing acute severity in the COVID disease. In the left-hand panel of To help reveal these interactions, we can colour by another feature. If we pass the whole explanation tensor to the color argument, the scatter plot will pick the best feature to colour by. In this case, it picks Fever Sympt (symptoms associated with fever) because that highlights that the the average age per category in the cohort has less impact on acute COVID severity for categories with a high Fever Sympt value.

The values of interaction between locally explained variables are a generalisation of those to higher order interactions. Fast exact computation of pairwise interactions is implemented for tree models. This returns a matrix for every prediction, where the main effects are on the diagonal and the interaction effects are offdiagonal. These values often reveal interesting hidden relationships, such as how the increased risk of death peaks for inpatients with mild febrile symptoms at the age between 20 and 34 years (see

Persistent homology identifies novel key proteomic features involved in severity

Based on the previous clinical characterisation of COVID severity, we exploited the proteomic plasma information available for the remaining 371 individuals in the cohort. Unfortunately, the particular geometry of inpatient's proteomes as embedded onto lower dimensional spaces resulted highly sensitive to parameter setups considered in downstream analyses according to entropy-based severity. In such a scenario, we computed topological invariant structures instead (see Video 5). These invariants, the so-called simplicial complex, qualitatively analyse features that persist across multiple scales. Such invariants can be classified over days 0, 3, and 7 by obtaining their generators through persistent homology (see Fig S5A). This analysis led us to identify unique protein configurations (see Figs S5B andS6) within inpatient proteomes based on their connected components [START_REF] Xia | Persistent homology analysis of protein structure, flexibility, and folding[END_REF][START_REF] Aktas | Persistence homology of networks: Methods and applications[END_REF]. Thus, the whole universe of proteome embeddings could be enclosed in the quotient space P=~under the given equivalent relation: for any p • 2P, p i ~pj if the projected proteome i "is similar to" on the set of all rotated tails, the so-called special orthogonal group SOð3Þ [START_REF] Hall | Lie Groups, Lie Algebras, and Representations[END_REF]. Hence, a partition with two classes of equivalence come out whose disjoint union determines the all three revealed groups of patients (see Video 6). Indeed, these classes of equivalence have as represents ½c b : = fx2P : x~c b g if b is a ball-like shape in SOð3Þ and ½c s : = fx2P : x~c s g if s is a startlike shape in SOð3Þ. These two classes can be visualised in upper and lower panels of Fig 4A andB. Now, to identify proteins whose profiles are invariants of each inpatient over the days of their stay f0; 3; and 7g, we primary analysed the classes of equivalence ½c b and ½c s of each partition by persistent homology per group. Specifically, we used persistence diagrams wherein we quantified the number of homology generators (see upper Fig 4C ) while testing their quality by means of confidence band generated by probabilistic boosting and density diffusion (see lower panel of Fig 4C). This enabled the identification of a unique set of proteins (see Video 1, Video 2, and Video 3) that encapsulated two-and threedimensional structures (Fig 4D andE) important to topologically characterise severity classification per group over the days of stay of each individual in the hospital. Thus, we mapped the transmembrane serine proteases TMPRSS5 and SS15 [START_REF] Huttlin | Architecture of the human interactome defines protein communities and disease networks[END_REF] as novel receptors taking part of the infection machinery. These proteins belong to the same family of TMPRSS2, a known receptor used by Sars-CoV-2, to enter the cell [START_REF] Hoffmann | SARS-CoV-2 cell entry depends on ACE2 and TMPRSS2 and is blocked by a clinically proven protease inhibitor[END_REF]. Both of those proteins were located amongst the dysregulated interactome of inpatients stratified as acute and strikingly also as mild. Furthermore, amid the proteins we found (see Video 1, Video 2, and Video 3 for the entire list), there were proteins in acute -BRK1, LAP3, SLC27A4, SLC39A14-, intermediate -SLC27A4-, and mild -BRK1, SLC27A4, and SLC39A14-inpatients functionally linked (Huang et al, 2009a(Huang et al, , 2009b) ) with AP3B1, BRD4, BRD2, CWC27, SLC44A2, and ZC3H18 whose profiles are reported to be likely involved in early infection caused by the virus [START_REF] Gordon | A SARS-CoV-2 protein interaction map reveals targets for drug repurposing[END_REF]. Finally, the functional analysis of the novel proteomic features resulted from our persistence analysis showed along with their ancestors two sharp clusters bound to pneumonia and inflammation pathways (see Fig S7).

Dynamic tracking of protein interactions required by the virus to efficiently infect the cell

Once we put the spotlight on individual proteins topologically important to discriminate COVID patients over time, we envisaged to capture their dynamics of functional interactions at regulatory levels. To this end, integrated protein-protein interaction networks were firstly constructed per group of patients using colocalization, coexpression, physical interactions, and shared domains [START_REF] Morilla | Assessment of protein domain fusions in human protein interaction networks prediction: Application to the human kinetochore model[END_REF][START_REF] Donaldson | The GeneMANIA prediction server: Biological network integration for gene prioritization and predicting gene function[END_REF]. Then, we enquired these graphs about their connection quality by means of degree and centrality distributions as shown in Fig 5 . Surprisingly, we could confirm (see Fig 5) that it neither was highly connected nor played an important modular role in the graphs. Next, to predict disease progression, we monitored the behavioural regulation of the nodes' graph aggregation on semi-supervised learning on a community composed by known and unknown protein interactions with ACE2 and TMPRSS2 [START_REF] Morilla | Deep models of integrated multiscale molecular data decipher the endothelial cell response to ionizing radiation[END_REF] (see the Video 1, Video 2, and Video 3, and Video 4). To compute such a tracking, we endowed the graphs with a tailored hybrid design (see Fig S8) covered in convolutional layers along with a spectral rule [START_REF] Defferrard | Convolutional Neural Networks on Graphs with Fast Localized Spectral Filtering[END_REF] as occurs in GCNs. The first model's performance yielded an accuracy, for each group of inpatients, of 0.49, 0.41, and 0.85 supported by 76, 146, and 355 samples regarding covid and noncovid feature representations, respectively (see Figs S9-S11 and Table 1). Those values increased to 0.71, 0.84, and 0.95, respectively, to the second learning model. We also found that their corresponding performances asymptotically tended to 0.65, 0.7, and 0.81 when the layers were largely increased from the 32 units in the convolutional architecture.

In that way, we learnt how ACE2 and TMPRSS2 interacted with the persistent novel candidates to explain the virus machinery at its entrance into the cell to put patients into mild, intermediate, or acute groups of severity over time (see Video 7, Video 8, Video 9, and Video 10). Hence, a primary set of proteins that led to acute severity consisted of the progressive aggregation of BCAN, CA2, CA12, CLEC4, FOLR1, FOLR2, IFNGR2, IGSF3 (R), ILR13A1(R), LAIR1, LRRN1, PCDH17, RTBDN, SEZ6L, SIGLEC6 (Schulte-Schrepping et al, 2020), and TNFRSF21 with respect to covid feature representation. Especially, CLEC4 belongs to a protein family (i.e., the C-type lectin receptor) involved in regulating immune reactivity through platelet degranulation whose expression significantly decreased in COVID-19 and correlated with disease severity [START_REF] Overmyer | Large-scale multi-omic analysis of COVID-19 severity[END_REF]. The interactions occurring early on during the infection amongst AXL, CD58, DDR1, DLK1, FCGR3A, TNFRSF12A, UXS1, and XPNPEP2 set the noncovid latent feature. Herein, we spotted CD58, a nonclassical monocyte, such as CD274 (PD-L1) known inhibitor of T-cell activation along with Arginase 1 [START_REF] Bronte | L-arginine metabolism in myeloid cells controls T-lymphocyte functions[END_REF][START_REF] Li | Arginine deficiency is involved in thrombocytopenia and immunosuppression in severe fever with thrombocytopenia syndrome[END_REF] highly expressed in neutrophils in COVID-19 patients or CD24 involved in neutrophil degranulation with an increased expression of neutrophil function [START_REF] Overmyer | Large-scale multi-omic analysis of COVID-19 severity[END_REF]. Therein, we also identified FCGR3A (encoding CD16a) that is regulating severitydependent alterations of the myeloid cell compartment during Sars-CoV-2 infection. Indeed, FCGR3A has been already found to be a nonclassical monocyte marker in COVID-19 (Schulte-Schrepping et al, 2020). Next, to the covid representations that determined intermediate severity of patients, we found the early aggregations of FR2, GALS4, IL1RN(R), ILR1(R), LRPAP1, RNF41, TRIM21(R), and VWV2. Remarkably, RNF41 plays a central role during interactions of Sars-CoV-2 with innate immune pathways because its interferon pathway is targeted by RNF41 (NSP15) [START_REF] Gordon | A SARS-CoV-2 protein interaction map reveals targets for drug repurposing[END_REF]. Then, the intermediate noncovid representations are governed by the interactions of CCR5 (intestinal pro-inflammatory), CPA1, LAMA2, PLA2G4A, PON3, SETMAR, TGFB1, and XCL1. Finally, aggregations of C4BPB, CD70(R), IPCEF1, MAVS(R), PLCG2, and THBS2(R) led to mild severity stratification to covid feature representations. At the same time, the interactions between CD200, MAPKAPK5 [START_REF] Kindrachuk | Antiviral potential of ERK/ MAPK and PI3K/AKT/mTOR signaling modulation for Middle East respiratory syndrome coronavirus infection as identified by temporal kinome analysis[END_REF], NTRK3(R), PRAP1(R), and XPNPEP2(R) set the noncovid feature representations. In these two lists, we might mention a similar effect on neutrophils and expression as CD58 to CD70 and CD200 [START_REF] Overmyer | Large-scale multi-omic analysis of COVID-19 severity[END_REF].

We checked that covid feature representation of patients with acute symptoms was functionally charaterised by a set of proteins involved in the fusion of virus enclosure to the host endosome membrane (GO:0039654) at the virus entrance into the cell. Overall, the interactions between ACE2 and TMPRSS2 and these persistent proteins were largely enriched in the immunoglobulin-like fold functional category. From a mere clinical stratification point of view, these conditions characterise most of the inpatients in the available clinical dataset suffering of cardiovascular complications. Noncovid representation of acute patients was strongly composed by membrane and signal peptide functional categories, protein tyrosine kinase and glycoprotein with extracellular and cytoplasmatic topological domains and transmembrane helix and integral component in the region biological processes. In this case, these conditions were felt on patients mainly suffering from diabetes of type 1 and immunosuppression.

Regarding the covid feature representation of intermediate patients, an overabundance of protein binding function is observed with diabetes of type 2 and normal variation diseases associated with such conditions. On the other side, the noncovid feature of this group is strongly enriched with disulfide bond and signal peptide functional categories. Therein, we found various terms directly linked with endosome viruses' machinery. Thus, we identified clathrin-dependent endocytosis (GO:0075512), host lysis, inhibition of host IKBKE, JAK1, RLR pathway, TBK1, and TLR pathway triggered by the virus in the host cell. Those symptoms were identified to a wide range of the complications described in the clinical dataset; in particular, asthma severity, chronic hepatitis C, immunosuppression after liver transplantation, diabetes, especially strong of type 2, heart and kidney complications, and hypertension.

Finally, the latent covid representation of mild patients' proteomes was functionally characterised by a weak overabundance of disulfide bond and glycosylation site (i.e., N-linked as GlcNAc, etc.). These categories were related with suppression by the virus of host-adaptative immune response (GO:0039504). Remarkably, The first set of figures before the "/" represents Model 1, whereas the second set of figures represents Model 2 with aggregation (see in Materials and Methods section: notes on the GCN). For the sake of simplicity, we only show acute subgroup patients. For further details, see Tables S1 andS2.

there were no disease-associated genes type-specific to these biological processes. The noncovid mild features were actually overrepresented by signal peptide, qualitatively similar to those features described to the intermediate noncovid patients. We will fully expose and discuss the intriguing implications of such results in the next section.

Discussion

Sars-CoV-2 has become in these two last years a real life-thread that has collapsed the health systems worldwide. Many efforts have been already done to structurally characterise the Sars-CoV-2 spike protein. To predict its severity, large mappings of proteins likely involved in the machinery applied by the virus to infect the cell have been reported [START_REF] Jackson | Mechanisms of SARS-CoV-2 entry into cells[END_REF][START_REF] Sokhansanj | Predicting COVID-19 disease severity from SARS-CoV-2 spike protein sequence by mixed effects machine learning[END_REF]. All these investigations have led to enormous advancements in COVID-19 treatment that consequently have given rise to efficient vaccines [START_REF] Kumari | A critical overview of current progress for COVID-19: Development of vaccines, antiviral drugs, and therapeutic antibodies[END_REF]. However, there is still some facets not well-characterised or yet sufficiently explored. In our attempt to contribute to this research, we computed an overall severity score based on WHO scales instrumental to provide a chart explaining the protein interactions required by the virus to stratify a patient's infection into mild, intermediate, or acute (Organisation world health, 2021). To this end, we made use of a double analysis linking symptoms to protein expressions and interactions with ACE2 and TMPRSS2. Merely from a stratification standpoint, we give novel information and verified known facts about COVID. As conclusion, we could claim that in itself, COVID-19 is not as harmful as it is in association with other risk factors such as age, febrile symptoms, or overweight. Indeed, most of the Covidpatients, though some were considered as outliers as indicated in the earlier sections, held a relative high entropy value of severity due to an eventual intubation, ventilation or supplementary oxygen requirement. To achieve those high peaks, Covid + should be overweight elderly people and present fever and/or respiratory symptoms. More importantly, we obtained functional evidence of how particular sequences of proteins interacted with the virus to block the immune systems during the infection. Thus, we could explain the infection fate towards the acute symptoms because an endosome acidification is produced during the infection initiating conformational proteins fusion [START_REF] Yang | SNX27 suppresses SARS-CoV-2 infection by inhibiting viral lysosome/ late endosome entry[END_REF]. In that way, Sars-CoV-2 could take advantage of such pathway to be endocytosed as happens with many types of viruses such as influenza A virus, alphaviruses, or HIV-1.

Regarding patients suffering from intermediate symptoms that endocytic process is caused over proteins that contain at least one coiled domain forming stiff bundles of fibres. Hence, proteins are modified by signaling upon creation of interchain disulfide bonds, which can produce stable, covalently linked protein complexes likely contributing to fold and stabilise proteins. In virus internalization, clathrin-mediated endocytosis could be then generated in response to getting assembled on the inside face of the cell membrane to cleave the host cell (CCV) by the action of DNM1/ dynamin-1 or DNM2/dynamin-2 (Jima & Hinshaw, 2018). Then, the virus may be delivering their content to early endosomes via CCV. These mechanisms could be expressing in Sars-CoV-2 using different ways as by lysing the host cell, blocking the host innate defenses via IKBKE/IKK-epsilon kinase inhibition, JAK1 protein, DDX58/RIG-I-like repector (RLR) what stabilises the antiviral state, TBK1 kinase inhibition to prevent IRF activation, or toll-like recognition receptor (TLR) pathway evasion, which makes the production of interferons to be inhibited and so to establish a stable antiviral state [START_REF] Chen | SARS-CoV-2 disrupts proximal elements in the JAK-STAT pathway[END_REF][START_REF] Aliyari | The evolutionary dance between innate host antiviral pathways and SARS-CoV-2[END_REF].

To mild cases, the Sars-CoV-2 protein could be preventing the tuned repertoire of self and nonself antigens' recognition of efficiently acting against the malicious effects of cell infection. In these cases, Sars-CoV-2 would be escaping the adaptive immune response by simple interference with the presentation of antigenic peptides at the surface of infected cells [START_REF] Sette | Adaptive immunity to SARS-CoV-2 and COVID-19[END_REF].

The use of plasma in COVID-19 research can be justified for many reasons, even if the cells of the immune system are not evaluated directly. Plasma contains a wide range of substances produced by the immune system, including antibodies and other proteins, which can provide valuable information about the body's response to the SARS-CoV-2 virus. By studying these substances in the plasma, researchers can gain insights into the immune system's response to the virus and how it may help or hinder the body's efforts to fight off the infection. Additionally, plasma can be collected and stored easily, making it a convenient and accessible source of samples for researchers studying COVID-19. Finally, plasma can be used in a wide range of experimental techniques, making it a versatile tool for researchers studying the immune system. However, it is important to note that the relationship between plasma findings and the immune system's response is not always clear, and further research may be needed to fully understand this relationship.

The approach presented in this work provides a more comprehensive and dynamic understanding of the molecular mechanisms of COVID-19 compared to more conventional methods such as antibody testing and PCR testing. These classic methods are limited in their ability to assess disease severity and provide a snapshot of the disease state, whereas the use of protein level measurements and clinical evaluation data, combined with GCNs, provides a dynamic understanding of the molecular interactions between proteins and can help determine the severity of a COVID-19 infection at an early stage.

Whereas traditional approaches rely on identifying the presence of the virus through antibodies or PCR, this study goes beyond that to provide insights into the molecular interactions and changes in protein expression that occur during the progression of the disease. This can aid in the identification of potential pharmacological treatments by studying the dynamic interactions between the most relevant proteins.

The use of GCNs in this study provides a unique and innovative approach to understanding the complex molecular mechanisms of COVID-19, offering valuable insights for future treatment and diagnostic tool development. However, this method may need further validation before becoming a routine diagnostic tool.

Overall, the results provided in this work contribute to gaining new insights into the mechanisms of the disease and how it affects the immune system and how nonlinear relationship between "message passing" proteins could particularly explain disease severity modulation during the Sars-CoV-2 infection.

Interpreting latent space of clinical features

The Shapley values are a widely used approach in cooperative game theory which have advantageous properties, as described by [START_REF] Yoshida | Verification of usefulness of SHAP values in interpretation of decision tree models[END_REF]. These properties facilitate a clear understanding of the computation and interpretation of Kapley values as explanations for machine learning models. This is achieved through an empirical approach using the shap Python package, which demonstrates the explanations for increasingly complex models [START_REF] Lundberg | A unified approach to interpreting model predictions[END_REF]. In this study, we apply Shapley values to explain the predictions made by a neural network model trained on clinical variables. The model was trained using the "adam" optimizer and "binary cross entropy" loss functions [START_REF] Chollet | GitHub[END_REF] to uncover why it makes different predictions for different individuals. To use SHAP, the model must take a 2D numpy array as input; so, a wrapper function was defined around the original Keras predict function [START_REF] Chollet | GitHub[END_REF]. To explain a single prediction (as shown in Fig 3A ), 50 samples from the dataset were selected to represent typical feature values and 500 perturbation samples were used to estimate the SHAP values, which required 500 * 50 evaluations of the model. To explain many predictions, the above process was repeated for 50 individuals. Note that as this explanation is based on a sampling approximation, each explanation may take a few seconds, depending on the machine setup (as shown in Fig 3B -D).

Notes on persistent homology (multi-topology)

Persistent homology is a mathematical theory used as a tool in topological data analysis (TDA) to study the topological features of a dataset (i.e., plasma proteome of each patient in the cohort). It is a type of algebraic topology that uses algebraic invariants to study the shapes of data. The main idea behind persistent homology is to study how the topological features of a dataset change as the scale of observation changes. This is often done by constructing a sequence of simplicial complexes from the data, each of which represents the data at a different scale, and then studying how the topological features of these simplicial complexes evolve over time.

In TDA, persistent homology is used to identify topological features of a dataset that are persistent, or stable, over a range of scales (i.e., multi-scale topology). These persistent features are considered to be the most significant topological features of the data, and they can be used to distinguish different classes of data or to identify patterns in the data. To use persistent homology in TDA, one typically starts by constructing a simplicial complex from the data, and then applying algebraic topology techniques to study the topological features of the complex. This can be done using specialised software tools, such as GUDHI, DIPHA, or Dionysus [START_REF] Maria | The gudhi library: Simplicial complexes and persistent homology[END_REF][START_REF] Gillani | DIONYSUS: Towards query-aware distributed processing of RDF graph streams[END_REF][START_REF] Tauzin | Giotto-tda: A topological data analysis toolkit for machine learning and data exploration[END_REF], which are designed to compute persistent homology and other topological invariants.

Notes on the GCN

We endowed the regulatory persistent-based protein networks with a convolutional design (i.e., neural networks) along with a spectral rule of node aggregation as in GCN [START_REF] Defferrard | Convolutional Neural Networks on Graphs with Fast Localized Spectral Filtering[END_REF]. The sequential combination of two hybrid models enabled the learning of interactions between ACE2 and TMPRSS2, and the persistent proteins needed by the virus to be spread in cells. Following this reasoning, we primarily made use of the identity matrix I as features and the adjacency matrix A (Model 1) contributing the model in the following spectral rule: where D is the degree matrix. Right after, we considered the metric given by distance of the shortest path to characterise the early aggregation of persistent proteins to ACE2 and TMPRSS2 as an additional feature in a second model (Model 2). Thus, we generated models with two layers computed by 32 units per layer and a 2D transformation of the activation function tanh. When applying the spectral rule, the relu activation function is applied at the beginning of the layer implementation instead of later on. The number of epochs was set to 250 and 5,000, respectively. We computed the stochastic gradient descent optimizer in the training task picking a learning rate and momentum regularization set to 0.001 and flagged true, respectively. The semi-supervised classification (covid or noncovid early proteins aggregation) of nodes in the persistent-based protein networks [START_REF] Kipf | Semi-Supervised Classification with Graph Convolutional Networks[END_REF] was performed by an in-house python script based on MXNet implementation (Faster Cheaper Leaner, 2021).

Figure 1 .

 1 Figure 1. Flowchart of the entire analysis. (A) Admission of inpatients presenting COVID symptoms and triage. (B) Blood extractions of inpatients at days 0, 3, and 7 after their admission. Measurement of 1,400+ plasma protein levels by proximity extension assay. (B1, B2, B3) Proximity extension assay protocol: immunoreaction, extension, and detection of levels by sequencing. (C) Record of clinical annotations and evaluation of discharged inpatients after 28 d of stay in the hospital. (D) The usage of algebraic invariants to study the shapes of each inpatient plasma proteome. (D1, D2, D3) Three stages of multi-topology analysis, namely, WHO scores of severity fitting and inpatient stratification based on that fit, protein candidates described by persistent homology per stratum, and prediction of disease progression based on convolutional neural networks constructed from those candidates and per stratum. Link https://www.olink.com/mgh-covid-study/. Modified and used with permission.

Figure 2 .

 2 Figure 2. HDBSCAN severity clustering by means of entropy measures. (A, B, C) UMAP projection of clinical data. (A) Stratification by Covid+ and Covid-. (B) Inpatients' stratification by WHO max score (i.e., a measurement correlated with the maximum WHO outcomes achieved by patients during their hospital stay and available in the clinical dataset). (C) Individual discrimination by Shannon's Entropy combined with Hdbscan clustering algorithm. The blank circles show inpatients considered as outliers by the dissimilarity of their symptoms. (D) Probability density function optimally fitted in accordance with Shannon's Entropy displaying three sharp peaks, namely, mild, intermediate, and acute associated with inpatients of the Covid cohort.

Figure 3 .

 3 Figure 3. Clinical evaluation of our entropy-based model on COVID severity score by the model and higher individual features. Initial explanations are based on a gradient boosted decision tree model trained on the COVID cohort. (A) Left: bar chart of the average SHAP value magnitude.Age was the most influential symptom, changing the predicted absolute COVID probability on average by two percentage points (0.02 on x-axis). Right: a set of beeswarm plots, where each dot corresponds to an inpatient in the cohort per significative symptom. The dot's position on the x-axis shows the impact that a symptom has on the model's prediction for a given inpatient. The piled-up dots mean the density of inpatients suffering from a symptom with similar impact on the model. Younger ages reduce the predicted Covid risk, elder ages increase the risk. (B) Globally stacked SHAP explanations clustered by explanation similarity. Inpatient profiles land on the x-axis. Red values increase the model prediction, blue ones decrease it. Two clusters stand out: On the left is a group with low predicted risk of suffering an acute Covid, whereas on the right, we have a group with a high predicted risk of suffering from acute COVID. (C) Top-bottom: locally stacked explanations clustered by explanation similarity for infection, lung, and respiratory symptoms. (D) Effect of a single feature across the whole cohort. Top-bottom: dependence plots for Age and Hypertension (HTN) features. These plots display inflection points in predicted age and hypertension as $Age cat$ and HTN (oldness by years on average and hypertension complaint per individual in the cohort) changes. Vertical dispersion at a single category of Age (resp. HTN) represents interaction effects with other features. To help reveal these interactions, we coloured by Fever (resp. BMI). We passed the whole explanation tensor to the colour argument in the dependence plots to pick the best feature to colour by. In this case, it selected fever symptoms (resp. Body Mass Index) because that highlights that the average age (hypertension) per inpatient has more (less) impact on Covid severity for categories with a low (high) Fever (BMI) value.

  Fig 3A, we observe this same effect in stacked red bars for our multi-class output task. To understand how a single feature affects the output of the model, we plotted the local value of that feature versus the value of the feature for all the inpatients in the clinical dataset (see Fig 3B). Now, we can zoom in some of those effects individually as shown in Fig 3C for the lung and respiratory (the one quantified as lowest in its contribution to severity learning) symptoms. Because locally explained values represent a feature's responsibility for a change in the model output, the plot in Fig 3D represents the change in predicted COVID severity as Age cat (the average age per category in the cohort), or preexisting hypertension change. Vertical dispersion at a single value of Age cat represents interaction effects with other features.

  Fig 3D -upper panel-), or that non-preexisting hypertension has less impact on individuals with a high BMI cat value (see Fig 3D -lower panel-).

Figure 4 .

 4 Figure 4. Multi-scale topology analysis flowchart. (A) Upper: class of equivalence [c b ] determined by umap projection of mild inpatients upon rotation on SO(3). Lower: umap projection of an inpatient's soluble proteome. (B) Upper: class of equivalence [c t ], taking as example to show the mild inpatient 101. Lower: umap projection of that inpatient's soluble proteome. (C) Upper: topological feature extraction from diagram of persistence of patient 101 and its later calibration. Lower: application of density diffusion for separating noise from robust signals in the persistence diagram of impatient 101. (D) Spotted loops of proteins enclosing dimension 2 structures important to explain severity stratification over time of patient 101. (E) Spotted voids of proteins enclosing dimension 3 structures important to explain severity stratification over time of patient 101.

Figure 5 .

 5 Figure 5. Regulatory gene networks predicting disease progression regarding ACE2 and TMPRSS2. (A) Mild patients. (B) Intermediate patients. (C) Acute patients. (D) Group of proteins nonfunctionally enriched in acute patients. Highlighted in red ACE2 as one initial seem in the downstream analysis of protein interactions occurred post-infection. Yellow enhances those proteins (ERBB3, CD48, CCR5, FCRL6, and PLA2G10, among others) with a higher connectivity degree in the networks.

Variable

  mass index: 0 = <18.5 (underweight) 1 = 18.5-24.9 (normal) 2 = 25.0-29.9 (overweight) 3 = 30.0-39.9 (obese) 4 = ≥40 (severely obese) 5 = Unknown HEART Pre-existing heart disease? HEART -(coronary artery disease, congestive heart failure, valvular disease) 0 = No 1 = Yes LUNG Pre-existing lung disease? LUNG -(asthma, COPD, requiring home O2, any chronic lung condition) 0 = No 1 = Yes KIDNEY Pre-existing kidney disease? KIDNEY -(chronic kidney disease, baseline creatinine >1.5, ESRD) 0 = No 1 = Yes DIABETES Pre-existing diabetes? DIABETES -(pre-diabetes, insulin and non-insulin dependent diabetes) 0 = No 1 = Yes HTN Pre-existing hypertension -HTN 0 = No 1 = Yes IMMUNO Pre-existing immunocompromised condition ? IMMUNO (active cancer, chemotherapy, transplant, immunosuppressant agents, aspenic) 0 = No 1 = Yes Resp_Symp Respiratory symptoms? Symp_Resp (sore throat, congestion, productive or dry cough, shortness of breath or hypoxia, or chest pain) 0 = No 1 = Yes Fever_Sympt Febrile symptom GI_Symp Any GI related symptoms at presentation (abdominal pain, nausea, vomiting, diarrhea) WHO 0 WHO score for day 0 study window -enrollment plus 24 h 1 = Death within 28 d 2 = Intubated, ventilated, survived to 28 d 3 = Non-invasive ventilation or high-flow nasal cannula 4 = Hospitalized, supplementary O2 required 5 = Hospitalized, no supplementary O2 required 6 = Not hospitalized WHO 3 WHO score for day 3 study window 1 = Death within 28 d 2 = Intubated, ventilated, survived to 28 d 3 = Non-invasive ventilation or high-flow nasal cannula 4 = Hospitalized, supplementary O2 required 5 = Hospitalized, no supplementary O2 required 6 = Not hospitalized WHO 7 WHO score for day 7 study window 1 = Death within 28 d 2 = Intubated, ventilated, survived to 28 d 3 = Non-invasive ventilation or high-flow nasal cannula 4 = Hospitalized, supplementary O2 required 5 = Hospitalized, no supplementary O2 required 6 = Not hospitalized WHO 28 WHO score on study day 28 1 = Death within 28 d 2 = Intubated, ventilated, survived to 28 d 3 = Non-invasive ventilation or high-flow nasal cannula 4 = Hospitalized, supplementary O2 required 5 = Hospitalized, no supplementary O2 required 6 = Not hospitalized WHO max WHO max category at 28 d (maximum WHO score within first 28 d with death being the maximum possible) 1 = Death within 28 d 2 = Intubated, ventilated, survived to 28 d 3 = Non-invasive ventilation or high-flow nasal cannula 4 = Hospitalized, supplementary O2 required 5 = Hospitalized, no supplementary O2 required 6 = Not hospitalized abs_neut_0_cat Absolute neutrophil count day 0 category: 1 = 0-0.99 2 = 1.0-3.99 3 = 4.0-7.99 4 = 8.0-11.99 5 = 12+ abs_lymph_0_cat Absolute lymphocyte count day 0 category: 1 = 0-0.49 2 = 0.50-0.99 3 = 1.00-1.49 4 = 1.50-1.99 5 = 2+ abs_mono_0_cat Absolute monocyte day 0 category 1 = 0-0.24 2 = 0.25-0.49 3 = 0.50-0.74 4 = 0.75-0.99 5 = 1.0+ creat_0_cat Creatinine day 0 category 1 = 0-0.79 2 = 0.80-1.19 3 = 1.20-1.79 4 = 1.80-2.99 5 day 3 category: 1 = 0-200 2 = 200-299 3 = 300-399 4 = 400-499 5 = 500+ abs_neut_7_cat Absulute neutrophil count day 7 category: 1 = 0-0.99 2 = 1.0-3.99 3 = 4.0-7.99 4 = 8.0-11.99 5 = 12+ abs_lymph_7_cat Absolute lymphocyte count day 7 category: 1 = 0-0.49 2 = 0.50-0.99 3 = 1.00-1.49 4 = 1.50-1.99 5 = 2+ (Continued on following page)

  S r ðA; IÞ i = å

Table 1 .

 1 Performance of the GCN hybrid model to mild COVID patients.

		Precision	Recall	f1-score	Support
	False	0.07/0.33	0.17/0.06	0.10/0.10	18
	True	0.95/0.95	0.89/0.99	0.92/0.97	337
	Accuracy			0.85/0.95	355
	Macro avg.	0.51/0.64	0.53/0.52	0.51/0.53	355
	Weighted avg.	0.91/0.92	0.85/0.95	0.88/0.93	355

Table 2 .

 2 Variable descriptions of the clinical dataset.

  Materials and MethodsSamplesData were provided by the Massachusetts General Hospital Emergency Department COVID-19 with Olink proteomic (publicly available at https://www.olink.com/mgh-covid-study/). Clinical dataset and plasma proteomes of 384 adult patients were distributed in 306 ð80%Þ patients that tested positive in COVID-19 that were named as Covid + and 78 patients tagged as Covidthat tested negative, even suffering from respiratory symptoms. Blood samples were drawn from Covid + patients on days 0, 3, and 7 after their admissions. The blood samples from Covidpatients were only drawn on day 0. The clinical dataset reported the 28 d outcome stratification according with the WHO scores of maximal severity described. Additionally, other clinical parameters including comorbidities such as patient age or preexisting medical complications (i.e., chronic kidney disease, any chronic lung condition such as asthma, diabetes) were also recorded. To anonymise this dataset, continuous variables were turned into categorical (see Table2and Supplemental Datas 1 and Supplemental Datas 2).

Table 2

 2

	. Continued
	Variable (subject_id)	Description (subject ID)

Proteome dynamics of COVID-19 severity Gauthier et al. https://doi.org/10.26508/lsa.202201624 vol 6 | no 5 | e202201624

Acknowledgements

We would like to thank the funding from the National Research Association (ANR) (Inflamex renewal 10-LABX-0017 to I Morilla), Consejería de Universidades, Ciencias y Desarrollo, fondos FEDER de la Junta de Andalucía (ProyExec_0499 to I Morilla), DHU FIRE Emergence 4, and the l'Agence de la Biom édecine.

Data Availability

All data used in this work are included in the article and/or supplementary material.

Ethics Statement

The Partners Human Research Committee approved the collection and analysis of samples and waived the requirement for informed consent.

Supplementary Information

Supplementary Information is available at https://doi.org/10.26508/lsa. 202201624.

Author Contributions

S Gauthier: data curation and formal analysis. A Tran-Dinh: methodology, and writing-original draft, review, and editing. I Morilla: conceptualization, data curation, formal analysis, supervision, funding acquisition, investigation, visualization, methodology, project administration, and writing-original draft, review, and editing.

Conflict of Interest Statement

The authors declare that they have no conflict of interest.