STAR★Methods

Key resources table

REAGENT or RESOURCE | SOURCE | IDENTIFIER
Chemicals, peptides, and recombinant proteins

Hexaphenoxycyclotriphosphazene (HPCTP) | Shanghai Macklin Biochemical Technology Co., Ltd. (Shanghai, China) | CAS: 1184-10-7
DMSO (cell culture grade) | Dalian Meilun Biotechnology Co., Ltd. (Dalian, China) | Cat.No.: PWL064
PBS buffer (1×) | Wuhan Pricella Biotechnology Co., Ltd. (Wuhan, China). | Cat.No.: PB180327
0.25% Trypsin solution | Wuhan Pricella Biotechnology Co., Ltd. (Wuhan, China). | Cat.No.: PB180225
Click Reaction Buffer | Shanghai Beyotime Biotechnology Co., Ltd. (Shanghai, China) | Cat.No.: C0071S

Critical commercial assays

Glucose uptake assay kit (2-NBDG fluorescence method) | Wuhan Servicebio Technology Co., Ltd. (Wuhan, China) | Cat.No.: G1744-200UL
Animal Tissue/Cell Total RNA Extraction Kit | Wuhan Servicebio Technology Co., Ltd. (Wuhan, China) | Cat.No.: G3640-50T
Cell RNA Extraction Kit | Absin Bioscience Inc. (Shanghai, China). | Cat.No.: abs60286-50T
Total RNA Extraction Reagent | Absin Bioscience Inc. (Shanghai, China). | Cat.No.: abs60154
SYBR Advanced qPCR SuperMix Kit | Absin Bioscience Inc. (Shanghai, China). | Cat.No.: abs60087-500T
One-step reverse transcription third-generation premix | Absin Bioscience Inc. (Shanghai, China). | Cat.No.: abs60246-100T
Antifade Mounting Medium | Coolaber Science & Technology Co., Ltd. (Beijing, China) | Cat.No.: SL1841

Deposited data

SEA (Similarity EnsembleApproach) | SEA Search Server | https://sea.bkslab.org/
TargetNet | TargetNet server | http://targetnet.scbdd.com
SwissTargetPrediction | Swiss Institute of Bioinformatics | http://www.swisstargetprediction.ch
ChEMBL database | EMBL-EBI | https://www.ebi.ac.uk/chembl/
Comparative ToxicogenomicsDatabase (CTD) | Comparative Toxicogenomics Database | https://ctdbase.org/
DAVID Bioinformatics Resources | Laboratory of Human Retrovirology and Immunoinformatics | https://davidbioinformatics.nih.gov/
Index | Index database | http://targetnet.scbdd.com/calcnet/index/
Protein DataBank (PDB) | RCSB Protein DataBank | https://www.pdbus.org/
PubChem | National Center for Biotechnology Information (NCBI) | https://pubchem.ncbi.nlm.nih.gov/

Experimental models: Cell lines

hCMEC/D3 complete culture medium | Zhejiang Meisen Cell Technology Co., Ltd. (Zhejiang, China) | Cat.No.: CTCC-003-0113-CM

Software and algorithms

ProTox 3.0 | Charité – Universitätsmedizin Berlin | https://tox.charite.de/protox3/
SWISS-MODEL | SWISS-MODEL | https://swissmodel.expasy.org/
ADMETlab 3.0 | Computational Biology and Drug Design Group, Xiangya School of Pharmaceutical Sciences | https://admetlab3.scbdd.com/
Maestro software | Schrödinger, LLC | Version 2024-4
Linear Constraint Solver (LINCS) algorithm. | GROMACS development team | GROMACS built-in algorithm
Particle-mesh Ewald (PME) method | GROMACS development team | GROMACS built-in algorithm
MM - GBSA method | Computational method | N/A
gmx_MMPBSA | gmx_MMPBSA developers | https://github.com/Valdes-Tresanco-MS/gmx_MMPBSA
GROMACS 5.1.4 software | GROMACS development team | Version 5.1.4
AutoDock 4.2.6 | The Scripps Research Institute | Version 4.2.6
AutoDockTools 1.5.6. | The Scripps Research Institute | Version 1.5.6

Other

PCR plates | Wuhan Servicebio Technology Co., Ltd. (Wuhan, China) | Cat.No.: PCR-9601W-HS

Experimental model and study participant details

Cell culture model

Human cerebral microvascular endothelial cells (hCMEC/D3) and the corresponding complete culture medium were obtained from Zhejiang Meisen Cell Technology Co., Ltd. (Zhejiang, China). Cells were cultured under standard conditions (37°C, 5% CO2) and used as an in vitro BBB model.

Chemical exposure model

Hexaphenoxycyclotriphosphazene was purchased from Shanghai Macklin Biochemical Technology Co., Ltd. (Shanghai, China).

Computational model

We obtained SMILES notations and two-dimensional structures of the six CTP compounds from the PubChem database. Among these, PFPCTP, HCCTP, HFCTP, EPFCTP, and HMCTP had directly available three-dimensional structures. For HPCTP, which lacked a 3D structure in the database, we generated an initial 3D model from its 2D structure using Avogadro-1.2.0 software and performed geometry optimization.

Method details

Establishment of criteria for assessing the potential of exogenous compounds to cross the blood - Brain barrier

Research has established that the passive diffusion of exogenous compounds across the blood-brain barrier (BBB) correlates strongly with specific physicochemical properties. Several early studies reported that physicochemical properties of compounds, including molecular weight, lipophilicity, hydrogen-bonding capacity, and polar surface area, were closely associated with their ability to cross the BBB via passive diffusion. These studies further explored the use of computationally predicted physicochemical parameters to preliminarily evaluate the passive BBB permeability of compounds. Meanwhile, the predicted permeability outcomes were initially compared with results obtained from in vivo animal experiments and in vitro cellular assays.26,27,42 Subsequent studies continuously expanded and refined the physicochemical descriptors associated with passive BBB diffusion. In addition, the validity and reliability of computational prediction approaches were further substantiated through comprehensive comparisons between predicted permeability results and experimental observations.25,28,38,43,44 Computational toxicology tools can reliably predict these properties currently. Available platforms include the QikProp module in Schrödinger software and ADMETlab 3.0.28,42,43,44

Latest literature lists the 12 computational physicochemical property parameters most relevant to the performance of compounds in passive diffusion through BBB, Based on this literature, we have built a scoring system to quantitatively assess BBB penetration possibility.44 This system incorporates 12 parameters obtained from Schrödinger software predictions. The selected parameters are: predicted BBB partition coefficient (QPlogBB), MDCK cell permeability, molecular weight, octanol/water partition coefficient (logP), dipole moment, molecular volume, hydrogen bond donor count, hydrogen bond acceptor count, Caco-2 cell permeability, predicted albumin binding affinity (logKhsa), polar surface area (PSA), and number of rotatable bonds. For each parameter, we identified literature-derived optimal ranges that favor passive BBB diffusion.25,26,27,28,29,44 Compounds earn one point for each parameter falling within its optimal range and the summation produces a composite BBB score ranging from 0 to 12 (Table 1).

Computer prediction of BBB permeability and toxicity for CTPs

Six CTP compounds molecular structures were converted to MOL format and imported into Maestro software (version 2024-4, Schrödinger) for subsequent processing. Structural processing employed the LigPrep module with the Optimized Potentials for Liquid Simulations (OPLS) force field. The optimized structures were then analyzed using the QikProp module, which generated 51 molecular descriptors for each compound. We selected 12 key parameters from this set for inclusion in our BBB penetration scoring system, as specified in Table 1.

The ADMETlab 3.0 platform, which implements machine learning and Quantitative Structure-Activity Relationship (QSAR) modeling, provided predictions for absorption, distribution, metabolism, excretion, and toxicity properties.30,43 We submitted SMILES notations of all six CTPs to this platform for analysis. From the resulting output parameters, we selected 10 descriptors for comparative analysis along with the 12 parameters obtained from QikProp predictions. The BBB index calculated by ADMETlab 3.0 served as an additional indicator for identifying compounds with enhanced BBB penetration potential.43,44,49 We next identified nomenclature inconsistencies between certain parameters shared by Schrödinger and ADMETlab 3.0 outputs. To resolve these discrepancies, we consulted official documentation from both platforms to establish accurate parameter correspondences and annotated these relationships accordingly. Additionally, we submitted the SMILES notations of all six CTPs into the ProTox 3.0 predictive toxicology platform. The resulting toxicity prediction radar charts provided supplementary data for our toxicity assessment.

Screening of potential target genes associated with CTPs and blood - Brain barrier function

Following initial screening, we employed five specialized databases to identify potential target genes for HPCTP, PFPCTP, HCCTP, and HFCTP. These databases included SEA, TargetNet, Index, SwissTargetPrediction, and ChEMBL. All searches were limited to Homo sapiens. We removed duplicate entries from the initial results and consolidated the remaining genes for each compound.

We then obtained high-confidence neurotoxicity-related genes from the Comparative Toxicogenomics Database (CTD). These genes were associated with neurological system disorders (NSD). The consolidated target genes from our database searches were cross-referenced with this neurotoxicity gene set. This intersection analysis identified overlapping target genes that may participate in both CTP exposure and neurological toxicity pathways.

Functional pathway analysis of target genes

The DAVID database provides an online bioinformatics resource for systematic functional annotation of gene and protein datasets. We employed this platform to analyze potential target genes related to blood-brain barrier (BBB) functions for the four selected CTPs. Our analysis included both Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment approaches. The GO examination encompassed three distinct categories: Biological Process (BP), Cellular Component (CC), and Molecular Function (MF). We established a statistical significance threshold of p-value <0.05 for identifying enriched functional terms. This analytical strategy helped clarify the biological significance and potential molecular mechanisms of the candidate genes.

Selection of core targets and molecular docking

Molecular docking (MD) simulates molecular interactions between ligands and biological targets. This computational method examines binding interactions with macromolecular proteins at the blood-brain barrier (BBB). The technique finds broad application in biomedical research, particularly in drug discovery and development. Combining molecular docking with other computational methods enhances prediction reliability and enables investigation of potential mechanisms for CTP penetration across the BBB.38,50

We identified potential solute carriers that interact with HPCTP by combining GO and KEGG enrichment results with established literature on BBB solute carriers. Three-dimensional structures of two solute carrier proteins, SLC2A1 (PDB ID: 6THA) and SLC6A3 (PDB ID: 9EO4), were obtained from the Protein DataBank. We removed crystallographic water molecules, heteroatoms, and native ligands from these structures, and missing residues and atoms were modeled using SWISS-MODEL. Protein preparation involved AutoDockTools 1.5.6. For both SLC2A1 and SLC6A3, we added hydrogen atoms, merged nonpolar hydrogens, and assigned Gasteiger partial charges. The prepared receptor structures were then converted to PDBQT format.51 For the HPCTP ligand, we added hydrogen atoms at physiological pH (7.(4) using Avogadro software. Subsequent processing in AutoDockTools involved merging nonpolar hydrogens and assigning partial charges before final format conversion.

Docking simulations employed AutoDock 4.2.6 with the Genetic Algorithm (GA). We conducted 100 independent runs for each docking experiment while maintaining default parameters for other settings. The resulting complex conformations underwent cluster analysis. We selected the highest-scoring representative conformation from each cluster for detailed interaction pattern examination.51

Molecular dynamics simulations

Molecular dynamics simulations offer detailed understanding of molecular interactions during chemical transport across the blood-brain barrier. These simulations help predict permeability involving active transport mechanisms.38 We conducted molecular dynamics simulations to evaluate binding stability between HPCTP and two solute carriers, SLC2A1 and SLC6A3. The simulations used GROMACS 5.1.4 software with the Amber ff99SB force field.52,53

We then generated force field parameters for HPCTP using the general Amber force field through AmberTools20. The complex system underwent solvation in a periodic dodecahedral box with the TIP3P water model. A minimum distance of 1.2 nm was maintained between solute atoms and box boundaries. The system was supplemented with 150 mM NaCl to simulate physiological conditions and establish charge neutrality. Energy minimization employed the steepest descent method until the system reached energy convergence.

The pre-equilibration phase consisted of consecutive 100 ps simulations using NVT and NPT ensembles, and with the temperature maintained at 310 K. We applied harmonic constraints of 1000 kJ mol−1 nm−2 to solute heavy atoms during this phase. Production simulations ran for 100 ns with all bond lengths constrained by the Linear Constraint Solver (LINCS) algorithm. Long-range electrostatic interactions utilized the particle-mesh Ewald (PME) method, while van der Waals interactions employed the Verlet Method. System temperature remained at 310 K through a velocity-rescaling thermostat, and pressure was controlled at 1 atm used a Parrinello-Rahman barostat.54,55,56

Trajectory analysis incorporated several measurements. Specifically, we calculated root-mean-square deviation for receptor Cα atoms and ligand heavy atoms relative to initial configurations using gmx rms. The radius of gyration, for each complex was determined through gmx gyrate, and intermolecular atomic contacts within 0.6 nm were quantified using gmx mindist. For binding free energy evaluation, we extracted 1000 equally spaced conformations from the equilibrated trajectory segment spanning 10–100 ns. Binding free energies were calculated with the MM - GBSA method in gmx_MMPBSA using default parameters.57

Computational Identification of Key HPCTP-Binding Residues in SLC2A1 and SLC6A3

Based on the per-residue energy decomposition analysis, the five residues of SLC2A1 that contributed most to HPCTP binding were individually mutated to alanine, yielding the T137A, P141A, H160A, Q161A, and W388A mutants. Similarly, the five most contributing residues of SLC6A3 were mutated to alanine to generate the R85A, L89A, T316A, H477A, and R544A mutants. For each mutant complex, 100 ns molecular dynamics simulations were performed, followed by binding free energy calculations, to evaluate the impact of these key residue mutations on the stability and binding affinity of SLC2A1 and SLC6A3 toward HPCTP.

Glucose uptake assay

The hCMEC/D3 cells were seeded into 12-well plates with 4×105 cells per well and cultured in complete medium at 37°C, 5% CO2 incubator for 24 h before exposure. HPCTP exposure concentrations were set at 0, 1, 10, and 100 nM. The 0 nM group contained 0.1% DMSO as the vehicle control.

After 24 h of exposure, 2-NBDG dye was added according to the manufacturer’s instructions. 2-2-NBDG is a fluorescently labeled glucose analog that is taken up by cells via glucose transporters (GLUTs). Following co-incubation of 2-NBDG with the cells for 6 h, fluorescent images were captured using a Carl Zeiss AG fluorescence microscope. Integrated Density values were quantified using ImageJ software. The results were subsequently normalized to obtain the relative fluorescence intensity ratios between the treatment groups and the control group.

Glucose permeability assay

6.5 mm Transwell chambers with 0.4 μm pore polycarbonate membranes were placed into 24-well plates, and hCMEC/D3 cells were seeded using the same procedure described above at a density of 2 × 104 cells per well, culturing for 48 h until cells complete confluence was achieved. Cells were exposed to HPCTP at concentrations of 0, 1, 10, and 100 nM for 24h with 400 μL exposure medium in the upper and lower layers of each Transwell insert. After exposure, 400 μL of glucose-containing complete medium (60 mM) was added to the upper chamber, while an equal volume of glucose-free complete medium was added to the lower chamber. Following 6h of incubation, medium from both the upper and lower chambers as well as cell lysates were collected for analysis, respectively. Glucose concentrations in the samples were determined using an enzyme-linked immunosorbent assay (ELISA). The Optical Density (OD) values of all samples were measured at 450 nm. A standard linear regression curve was generated by plotting the concentrations of the standards against their corresponding OD values, and glucose concentrations in the samples were subsequently calculated based on this calibration curve. The obtained concentrations were normalized to generate relative percentage values for the treatment groups compared with the control group.

Reverse transcription quantitative polymerase chain reaction

Total RNA from hCMEC/D3 cells was extracted according to the instructions of the Servicebio Cell Total RNA Extraction Kit. Reverse transcription quantitative polymerase chain reaction (RT-qPCR) was performed according to the reagent kit instruction manual. Amplification was carried out using a Bio-Rad PCR system with an annealing temperature of 60 °C. The human-specific primers for SLC2A1 mRNA were as follows: forward, 5′-TGAGCATCGTGGCCATCTTT-3′; reverse, 5′-CCGGAAGCGATCTCATCGAA-3′. GAPDH was used as the endogenous control, with the following primers: forward, 5′-AATGGGCAGCCGTTAGGAAA-3′; reverse, 5′-GCGCCCAATACGACCAAATC-3′.

Quantification and statistical analysis

Statistical analyses were performed using SPSS. Comparisons among multiple groups were performed using one-way analysis of variance (ANOVA). Significance levels were denoted by ∗p < 0.05, and standard error was presented as mean ± standard error of the mean (SEM).