Investigating and Modelling an Engineered Millifluidic In Vitro Oocyte Maturation System Reproducing the Physiological Ovary Environment in the Sheep Model

In conventional assisted reproductive technologies (ARTs), oocytes are in vitro cultured in static conditions. Instead, dynamic systems could better mimic the physiological in vivo environment. In this study, a millifluidic in vitro oocyte maturation (mIVM) system, in a transparent bioreactor integrated with 3D printed supports, was investigated and modeled thanks to computational fluid dynamic (CFD) and oxygen convection-reaction-diffusion (CRD) models. Cumulus-oocyte complexes (COCs) from slaughtered lambs were cultured for 24 h under static (controls) or dynamic IVM in absence (native) or presence of 3D-printed devices with different shapes and assembly modes, with/without alginate filling. Nuclear chromatin configuration, mitochondria distribution patterns, and activity of in vitro matured oocytes were assessed. The native dynamic mIVM significantly reduced the maturation rate compared to the static group (p < 0.001) and metaphase II (MII) oocytes showed impaired mitochondria distribution (p < 0.05) and activity (p < 0.001). When COCs were included in a combination of concave+ring support, particularly with alginate filling, oocyte maturation and mitochondria pattern were preserved, and bioenergetic/oxidative status was improved (p < 0.05) compared to controls. Results were supported by computational models demonstrating that, in mIVM in biocompatible inserts, COCs were protected from shear stresses while ensuring physiological oxygen diffusion replicating the one occurring in vivo from capillaries.

Introduction
In vitro maturation (IVM) of oocytes is an assisted reproductive technology (ART), and an alternative to controlled ovarian hyperstimulation (COH). In COH, patients undergo pituitary downregulation followed by exogenous gonadotropin administration and retrieval of in vivo matured oocytes [1]. However, in IVM, meiotically immature oocytes are retrieved from unstimulated ovaries and cultured in vitro prior to fertilization under the influence of different additives, including but not limited to gonadotrophins. IVM has great potential to become the technique of choice in human and animal reproductive medicine, even though, despite its long history, it needs technical improvements.
In animal ARTs, it is the predominant method of maturation of oocytes used for in vitro production of embryos for commercial and research purposes [2]. The first report of oocyte isolation from ovarian follicles and IVM was published by Pincus and Enzmann

Collection of Ovaries and COC Retrieval
Ovaries from prepubertal lambs (less than 6 months of age) were recovered at local slaughterhouses, from the animals subjected to routine veterinary inspection in accordance with the specific health requirements stated in Council Directive 89/556/ECC and subsequent modifications. Ovaries were transported to the laboratory at room temperature within 4 h after collecting. For COC retrieval, ovaries underwent the slicing procedure [57]. Follicular contents were released in sterile Petri dishes containing phosphate buffered saline (PBS) and observed under a Nikon SMZ18 stereomicroscope equipped with a transparent heating stage set up at 37 • C (Okolab S.r.l., Napoli, Italy). Only COCs with at least three intact cumulus cell layers and homogenous cytoplasm were selected for culture [56].

Native Millifluidic IVM (Native mIVM)
In the first experiment, mIVM was performed by inserting COCs as they are (in the further text "native COCs"), on the bottom of the culture chamber of a commercial LB1 bioreactor. This culture system has been called "native mIVM" in the further text ( Table 1). The LB1, consisting of a polydimethylsiloxane (PDMS) chamber with a transparent top and bottom, was designed to reproduce the typical volume of the single well of a 24-well plate (2 mL) with a flow inlet and an outlet for the perfusion of culture media, was assembled as described by Mastrorocco et al. (2021) [56]. Briefly, COCs were placed inside the cell chamber after adding 300 µL of IVM medium. Then, the chamber was hermetically sealed and filled with an additional 1700 µL of IVM medium through the inlet tube. The bioreactor was kept for two hours in the incubator to allow the COCs to stabilize at the bottom of the chamber before flow activation. Two hours after COC loading, the culture chamber was connected to a millifluidic circuit composed of a mixing chamber, a reservoir of 18 mL with IVM medium, and a peristaltic pump which permits the circulation of the culture medium. A silicone tube connected with a 0.22 µm filter, in the mixing chamber, allowed medium oxygenation. The flow rate was temporally set at 450 µL/min for the time to quickly fill the circuit and then regulated to 50 or 100 µL/min, according to the experimental design. Flow rates were chosen based on previous studies in the sheep which reported that ovarian blood flow gradually decreases from the order of milliliters/min in the ovarian artery [51,52], to microliters/min in the capillary network around the developing follicles, varying in relation to various aspects such as the anatomical district and the ovarian cycle stage [53,54]. Table 1. Nomenclature of used millifluidic IVM (mIVM) culture systems.

Control
Static IVM in 4-well dish with no inserts.

Native
Dynamic mIVM in bioreactor without inserts.

Ring
Dynamic mIVM in bioreactor with and ring support placed at the bottom of the bioreactor chamber.

Concave + ring
Dynamic mIVM in bioreactor with concave support placed at the bottom of the bioreactor chamber and the ring support placed on its top.

Concave + ring + plain = Box
Dynamic mIVM in bioreactor with a box obtained by sequentially assembling the concave support at the bottom of the bioreactor chamber, the ring on its top and the plain support used as a lid.
Concave + ring + alginate = Alginate-filled supports Dynamic mIVM in bioreactor with COCs embedded in alginate gel loaded in concave support and ring inserted at the bottom of the bioreactor chamber.

mIVM with 3D Printed Biocompatible Supports
In the second experiment, mIVM culture was performed by inserting COCs within biocompatible supports placed at the bottom of the bioreactor chamber in order to recreate a protecting structure around them, similar to that of an ovarian follicle or part of it ( Table 1). The supports, consisting of a ring (inner/outer radius = 5.5/7 mm) and a plain or a concave disk (outer radius 7 mm, maximum depth = 1.5 mm), were designed with Fusion 360 and fabricated at the Research Center 'E. Piaggio', University of Pisa, using a stereolithographic 3D printer (Form2, FormLabs) loaded with a biocompatible photopolymeric material (Dental SG resin, Formlabs) [59]. The configurations tested were: (i) ring; (ii) concave + ring; (iii) concave + ring + plain (namely ring, concave + ring and "box" in the further text; Table 1; Figure 1). In all conditions, COCs were uploaded in the chamber after inserting the supports and adding the first 300 µL of IVM medium. Next, the chamber was then completely filled, closed, and connected to the pump as described above. The flow rate was set at 50 µL/min.

mIVM in Biocompatible Supports Filled with Alginate Gel
In the third experiment, a concave support and a ring were inserted at the bottom of the bioreactor chamber and filled with alginate gel with the aim to reproduce 3D conditions for COC culture. Inside the support, 150 µL of 1% w/v sodium alginate pre-equilibrated in IVM medium was released [55] and, within it and appropriately spaced, 20-25 COCs were uploaded. This culture condition was named "alginate-filled support" (Table 1). A 100 mM calcium chloride (CaCl 2 ) solution was then sprayed onto the surface of the gel at the top edge of the concave support [60]. After this step, CaCl 2 was added in droplets to allow complete alginate jellification and COC trapping. After 5 min, CaCl 2 was removed. The bioreactor was filled with IVM medium, sealed and kept for two hours in the incubator. After that, it was connected to the reservoir of IVM medium, and a pump set to a flow rate of 50 µL/min started the millifluidic culture as described above. In each replicate, 20-25 COCs were cultured. After 24 h of IVM culture, sodium alginate was removed by calcium chelation with 2% w/v sodium citrate in IVM medium for 5 min at 38.5 • C under 5% CO 2 in air.

mIVM in Biocompatible Supports Filled with Alginate Gel
In the third experiment, a concave support and a ring were inserted at the bottom of the bioreactor chamber and filled with alginate gel with the aim to reproduce 3D conditions for COC culture. Inside the support, 150 μL of 1% w/v sodium alginate preequilibrated in IVM medium was released [55] and, within it and appropriately spaced, 20-25 COCs were uploaded. This culture condition was named "alginate-filled support" (Table 1). A 100 mM calcium chloride (CaCl2) solution was then sprayed onto the surface of the gel at the top edge of the concave support [60]. After this step, CaCl2 was added in droplets to allow complete alginate jellification and COC trapping. After 5 min, CaCl2 was removed. The bioreactor was filled with IVM medium, sealed and kept for two hours in the incubator. After that, it was connected to the reservoir of IVM medium, and a pump set to a flow rate of 50 μL/min started the millifluidic culture as described above. In each replicate, 20-25 COCs were cultured. After 24 h of IVM culture, sodium alginate was removed by calcium chelation with 2% w/v sodium citrate in IVM medium for 5 min at 38.5 °C under 5% CO2 in air.

Assessment of Cumulus Expansion and Oocyte Denuding
After all three types of experiments, COCs were recovered and cumulus expansion was checked. COCs showing cumuli with continuous edges, consisting of cells in close contact each other, were classified as compact, whereas cumuli showing discontinuous edges following cell detachment and production of a viscous extracellular matrix were classified as expanded. The percentage of expanded COCs was noticed, as it represents a response of immature oocytes to the presence of gonadotropins in the medium even if it does not fully ensure that the maturation is achieved. COCs underwent cumulus cell removal by incubation in TCM-199 with 20% FCS containing 80 IU hyaluronidase/mL and aspiration in and out of finely drawn glass pipettes. Denuded oocytes were evaluated for meiotic stage and matured ones were used to assess bioenergetic/oxidative status.

Oocyte Mitochondria and ROS Staining
Oocytes were washed three times in PBS with 3% BSA and incubated for 30 min in the same medium containing 280 nmol/L MitoTracker Orange CMTM Ros (Thermo Fisher Scientific, Waltham, MA, USA) at 38.5 • C under 5% CO 2 in air. After incubation with MitoTracker, oocytes were washed in PBS with 0.3% BSA and incubated for 15 min, at 38.5 • C under 5% CO 2 , in air in the same medium containing 10 µmol/L 2 ,7dichlorodihydrofluorescein diacetate (H 2 DCF-DA) to detect the dichlorofluorescein (DCF) and localize intracellular sources of ROS [61]. After incubation, oocytes were washed in PBS without BSA and fixed overnight at 4 • C in 4% paraformaldehyde (PFA) solution in PBS [62]. Particular attention was applied to avoid sample exposure to the light during staining and fixing procedures and to reduce photobleaching.

Oocyte Nuclear Chromatin Evaluation
To evaluate oocyte nuclear chromatin, after the fixation in 4% PFA in PBS, oocytes were stained with 2.5 µg/mL Hoechst 33,258 in 3:1 (v/v) glycerol/PBS mounted on microscope slides with coverslips, sealed with nail polish, and kept at 4 • C in the dark until observation. Slides were examined under the epifluorescence microscope (Nikon Eclipse 600; Nikon Instruments, Firenze; ×400 magnification) equipped with a B-2A (346 nm excitation/460 nm emission) filter. Oocytes were evaluated in relation to their meiotic stage, and classified as germinal vesicle (GV), metaphase to telophase I (MI to TI) and MII with the first polar body extruded [63]. Oocytes showing either multipolar meiotic spindle, irregular chromatin clumps, or absence of chromatin were considered as abnormal [55].

Assessment of Mitochondrial Distribution Pattern and Intracellular ROS Localization
Oocytes at the MII stage were observed at ×600 magnification in oil immersion with a Nikon C1/TE2000-U laser scanning confocal microscope (Nikon Instruments, Firenze). A 543 nm helium/neon laser and the G-2A filter were used to detect the MitoTracker Orange CMTM Ros (551 nm excitation and 576 nm emission). A 488 nm argon ion laser and the B-2A filter were used to detect DCF (495 nm excitation and 519 nm emission). Scanning was conducted with 25 optical sections from the top to the bottom of the oocytes, with a step size of 0.45 µm to allow 3D distribution analysis. The mitochondrial distribution pattern was evaluated on the basis of previous studies: (1) finely granular, with small mitochondria aggregates spread throughout the cytoplasm, typical of immature oocytes; (2) perinuclear and subplasmalemmal (P/S) distribution of mitochondria forming large granules, which is indicator of cytoplasmic maturity; (3) abnormal, with irregular distribution of mitochondria [62]. Concerning intracellular ROS localization, oocytes with intracellular ROS distributed throughout the cytoplasm, together with areas/sites of mitochondria/ROS overlapping, were considered healthy.

Quantification of Mitochondrial Activity, Intracellular ROS Levels, and Mitochondria-ROS Colocalization
In each individual oocyte, MitoTracker and DCF fluorescence intensities were measured at the equatorial plane and at the excitation/emission, as described above by use of EZ-C1 Gold Version 3.70 image analysis software platform for Nikon C1 confocal microscope. A circular area was drawn in order to measure only the region including cell cytoplasm. The fluorescence intensity within the programmed scan area was recorded and plotted against the conventional pixel unit scale (0-255). Mitochondrial activity and intracellular ROS levels were recorded as MitoTracker Orange CMTM Ros and DFC fluorescence intensity in arbitrary densitometric units (ADU). Parameters related to fluorescence intensity, such as laser energy, signal detection (gain), and pinhole size, were maintained at constant values for all measurements. The degree of mitochondria-ROS colocalization, reported as a biomarker of healthy oocytes [62,63] was quantified by the overlap coefficient between MitoTraker Orange CMTM Ros and DCF fluorescence intensity signals.

Computational Models of Millifluidic IVM Systems
To evaluate the effects of flow on COC in vitro cultures, computational models of tested culture systems were developed in collaboration with the Research Center 'E.
Piaggio' of the University of Pisa. The models were implemented in COMSOL Multiphysics (Stockholm, Sweden), coupling incompressible Navier-Stokes and CRD equations. The bioreactor was modeled as a cylindrical chamber with PDMS lateral walls (oxygen permeability = 3.78 × 10 −11 mol m m −2 s −1 mmHg) [64,65], while the static well was implemented as a cylindrical geometry with a media layer of 2.1 mm covered by an oil layer of 2.1 mm. The support geometries were imported from the Fusion 360 files and placed at the bottom of the bioreactor reproducing the different configurations investigated. The fluid viscosity of the culture medium was assumed equal to that of water at 37 • C (i.e., 0.6913 mPa × s), while the oxygen diffusion coefficient was assumed equal to 3 × 10 −9 in the media and 2.5 × 10 −9 m 2 /s in the alginate gels [66] and to 2 × 10 −9 m 2 /s in the oil layer [67]. Finally, COCs were modeled as spheres with 0.3 mm radius (number = 21) and a consumption rate (OCR = 4 pmol/min) estimated from data reported in the literature for bovine COCs [68]. The initial oxygen concentration and the oxygen partial pressure were set at 0.21 mol/m 3 and 159 mmHg, respectively.

Statistical Analysis
The proportions of oocytes showing different chromatin configurations and mitochondria distribution patterns were compared between groups by the Chi-square test. Mitochondria and ROS quantification analysis was conducted on oocytes at MII stage. Data (mean ± standard deviation [SD]) were compared by one-way ANOVA followed by Dunnett's Multiple Comparison Test. Differences with p < 0.05 were statistically significant.

Native mIVM Did Not Sustain Oocyte Nuclear and Cytoplasmic Maturation
In experiment 1, in three independent runs, 229 COCs were cultured and 222 of them were evaluated. Under both flow rates of 50 and 100 µL/min, after 24 h mIVM, COCs were found at the edges of the bioreactor chamber. The cumulus expansion rate was 25% compared to 100% of COCs cultured in control static IVM. In mIVM, the percentage of oocytes which reached the MII stage was significantly lower than in controls (p < 0.001, Table 2). Conversely, the percentage of oocytes that remained at the MI-TI stages was significantly increased either at 50 µL/min (p < 0.01, Table 2) or at 100 µL/min (p < 0.001, Table 2). Moreover, a statistically significant increase of the percentage of oocytes showing abnormal chromatin configurations, such as multipolar meiotic spindle or irregular chromatin clumps, was observed (p < 0.05, Table 2).  In most of the control MII oocytes, the mitochondria distribution was P/S, which is a biomarker of cytoplasmic maturity (Table 3). Instead, after mIVM, either under 50 or 100 µL/min flow rate, significantly higher proportions of MII oocytes with finely granular mitochondria distribution were found (p < 0.05 and p < 0.01, respectively; Table 3 and Figure 2) with corresponding reduction of the percentage of oocytes with healthy P/S pattern (p < 0.05 and p < 0.001, respectively; Table 3 and Figure 2). At both tested conditions, oocyte mitochondrial membrane potential, ROS levels, and mitochondria/ROS colocalization were significantly reduced compared with controls (p < 0.001; Figure 2).

mIVM with Open Biocompatible Supports Sustain Oocyte Nuclear and Cytoplasmic Maturation
Based on the results of experiment 1, in Experiment 2, the effects of including the COCs in biocompatible supports placed inside the culture chamber, were tested. The oocytes were cultured in six replicates. In each replicate, at least one or two test conditions were examined and compared to the controls. A total of 372 COCs were cultured and 348 COCs were analyzed. By using the ring support, after mIVM, cumulus expansion was 50% compared with 100% of control COCs. The percentage of oocytes reaching the MII stage was lower compared to controls (p < 0.05; Table 4). When the concave support and the ring support were tested, cumulus expansion reached 80% and the percentage of matured oocytes returned to be comparable to that of static controls. Instead, in the box support, cumulus expansion was reduced (20%) and oocyte nuclear maturation was inhibited (p < 0.001, Table 4). Figure 3 shows COCs cultured under native mIVM (A), with the ring (B) and the concave + ring (C) supports.  Both the mature oocytes obtained by mIVM culture in the presence of ring and concave + ring supports showed healthy P/S type mitochondrial distribution pattern comparable to the controls (Table 5). On the contrary, after mIVM in the box support, the percentage of matured oocytes with P/S mitochondrial distribution pattern was significantly reduced (p < 0.05, Table 5), with a significantly higher percentage of matured oocytes showing abnormal, particularly large mitochondria granulations spread all over the cytoplasm or located in specific cytoplasmic domains (p < 0.001, Table 5 and Figure 4). In oocytes matured with the ring support, mitochondrial activity, ROS levels, and mitochondria/ROS colocalization were significantly reduced compared to controls (p < 0.001; Figure 4). In MII oocytes derived from mIVM with concave + ring, no differences were noticed for mitochondrial activity and ROS levels compared with controls. However, as both values were higher than controls, a statistically significant increase in mitochondria/ROS colocalization was found (p < 0.001; Figure 4). For oocytes matured in the box support, quantification parameters did not differ with those of controls.  Both the mature oocytes obtained by mIVM culture in the presence of ring and concave + ring supports showed healthy P/S type mitochondrial distribution pattern comparable to the controls (Table 5). On the contrary, after mIVM in the box support, the percentage of matured oocytes with P/S mitochondrial distribution pattern was significantly reduced (p < 0.05, Table 5), with a significantly higher percentage of matured oocytes showing abnormal, particularly large mitochondria granulations spread all over the cytoplasm or located in specific cytoplasmic domains (p < 0.001, Table 5 and Figure 4).

D-mIVM in Alginate Layer Sustain Oocyte Nuclear and Cytoplasmic Maturation
In Experiment 3, the effects of mIVM with COCs included in supports filled with an alginate layer were analyzed. A total of 150 oocytes were processed in three independent runs and 146 of them were analyzed. After 24 h of IVM in alginate-filled supports, cumulus expansion ranged between 80 and 100%, and the percentage of MII oocytes was comparable to that of controls (Table 6). Additionally, no differences in mitochondria distribution pattern were found compared to controls (Table 7). These MII oocytes showed increased mitochondrial activity (p < 0.05), ROS levels (p < 0.001; Figure 5), and mitochondria/ROS colocalization (p < 0.05; Figure 5). Figure 6 shows representative photomicrographs of oocytes cultured in the analyzed condition. It can be seen that oocytes matured in alginatefilled supports showed P/S mitochondria distribution pattern and increased mitochondriaand ROS-related fluorescence intensities.
In the ring and concave + ring configurations, shear stress values were comparable with the bioreactor set at 50 μL/min in native conditions (3.14 × 10 −6 and 3.28 × 10 −6 Pa) and average oxygen concentrations were not significantly increased with respect to static conditions (0.184 and 0.182 mol/m 3 , respectively). Differently, in the box configuration, despite protection of COCs from flow and shear stress, oxygen concentration was lower (0.146 mol/m 3 ).
Chi square test: Not significant.
In the ring and concave + ring configurations, shear stress values were comparable with the bioreactor set at 50 µL/min in native conditions (3.14 × 10 −6 and 3.28 × 10 −6 Pa) and average oxygen concentrations were not significantly increased with respect to static conditions (0.184 and 0.182 mol/m 3 , respectively). Differently, in the box configuration, despite protection of COCs from flow and shear stress, oxygen concentration was lower (0.146 mol/m 3 ).
Finally, in the alginate-filled support, due to the presence of alginate, direct flow rate was absent and so the shear stress and average oxygen levels resulted as being equal to 0.182 mol/m 3 . Figure 7 shows the CFD and oxygen CRD models and relative data for all the examined conditions. Finally, in the alginate-filled support, due to the presence of alginate, direct flow rate was absent and so the shear stress and average oxygen levels resulted as being equal to 0.182 mol/m 3 . Figure 7 shows the CFD and oxygen CRD models and relative data for all the examined conditions.

Discussion
The first step of this study was to investigate whether the native mIVM, without any kind of protective or 3D-inclusive support, could be advantageous for COC maturation.

Discussion
The first step of this study was to investigate whether the native mIVM, without any kind of protective or 3D-inclusive support, could be advantageous for COC maturation. This hypothesis was formulated on the basis on previous studies in other cell systems, which showed that millifluidic cell cultures improve nutrient and oxygen supply [65,[69][70][71].
To the best of our knowledge, this is the first study on native mIVM in sheep, a large animal model. In a previous study from our unit, we found that mIVM, not even increasing the maturation rate, boosted oocyte bioenergetic/oxidative status compared with the static conditions. However, in that study, a 3D-IVM system was used by constructing COCincluding alginate microbeads via bioprinting technologies [56]. Differently, in this study, we wanted to test direct effects of mIVM itself on the COC. To our knowledge, this is the first study analyzing the effect of dynamic mIVM on immature oocytes of farm animals not submitted to superovulation treatments but derived from the spontaneous ovarian cycle, which are suitable models for fertility rescue or preservation programs in both humans and in animals. Previous studies were performed in oocytes from superovulated mice [47,48] or in bovine oocyte by using a microfluidic but not dynamic system [49]. The results of the first experiment clearly indicated that exposure of native COCs to the tested flow conditions did not support oocyte nuclear and cytoplasmic maturation, as most of the oocytes did not reach the MII stage, and a few matured oocytes showed compromised bioenergetic/oxidative status.
In order to evaluate the reasons for such significant oocyte damage, computational models of the tested culture systems were developed. As shown by the CFD models, the shear stress on COCs was found to increase correspondingly to the flow rate. Moreover, according to the CRD computation, higher average oxygen concentrations were found with respect to static well controls. Thus, these conditions may have negatively affected COC maturation, and may be responsible for observed inhibition of oocyte meiotic maturation and mitochondrial unbundling and reduced activity. Indeed, in vivo, COCs are not directly exposed to blood flow and nutrients, and oxygen reaches them by diffusion from capillaries. It is worth mentioning that ovaries have a complex blood supply system required to support ovarian functions. For example, in the sheep model as well as in humans, growing follicles are found in the corticomedullary border, a region particularly well-supplied with blood vessels [72]. Follicle microvasculature never passes beyond the basal lamina before ovulation when it invades into the granulosa layer [73,74]. The ovarian vasculature is involved in oxygen tension regulation as oxygen released from capillaries first reaches the thecal cells and then lower amounts diffuse through the basement membrane to reach the multiple layers of granulosa cells and the COC. The oxygen concentration at the follicle surface in developing human follicles is predicted with a range of 0.121-0.131 mol/m 3 [50,75]. Due to the lack of direct blood supply, as the diameter of ovarian follicles greatly increases during follicle maturation, the diffusion distance for gasses inside the follicle also increases, and this led to a continuous decrease of oxygen concentration in the follicular fluid during follicular maturation, reaching the lowest levels in preovulatory follicles [50,76,77]. These data are in line with compromised COC maturation found under flow rate and oxygen concentrations higher than those of physiological conditions. Based on the results of experiment 1 and information of in vivo condition reported in the literature, in experiment 2, we evaluated whether the insertion of the COCs in biocompatible supports, placed inside the culture chamber, could protect the COCs from the shear stress and from possible side effects of high oxygen levels, maintaining the advantages of millifluidics (i.e., continuous arrival of new nutrients, the removal of metabolism waste). By inserting the ring support, less significant cell damage was observed. In fact, the percentages of nuclear maturation, although reduced compared to the controls, were almost high as absolute value (53%) and significantly higher than that of IVM culture in the absence of supports at the same flow rate (53% in Table 4 vs. 27% in Table 2 at 50 µL/min; p < 0.05). Furthermore, the good P/S mitochondrial pattern was preserved in a percentage comparable to the controls, although the quantification data indicated a significant trend towards loss of viability. More evident improvement was obtained by inserting two supports, the concave one and the ring one, as the nuclear maturation and P/S mitochondrial pattern were reached at percentages comparable to the controls. Additionally, quantification data of mitochondria activity and ROS levels tended to increase, indicating cellular activity and viability, even if they did not attain statistical significance, which was observed for colocalization. Conversely, the box configuration significantly inhibited oocyte nuclear maturation and mitochondria distribution, although quantification parameters were not affected. Computational models revealed that, in the support configurations (ring and concave + ring), although shear stress values were comparable with those of the native mIVM culture at 50 µL/min without supports, oxygen concentration did not significantly increase. Moreover, in both configurations, at the end of the culture time COCs were found to be much closer together. These two conditions, lower oxygen levels and increased proximity between COCs, may have had a positive impact on COC culture conditions (e.g., by promoting the exchange of signaling molecules). In addition, the concave support may also have had a positive effect replicating in vivo like curvatures as in previous studies [78]. Differently, in the box configuration, despite COCs being protected from flow, oxygen concentrations were lower. Indeed, the box configuration had a low oxygen permeability which may have inhibited oxygen diffusion towards the COCs, thus affecting their major cellular activities. Moreover, the box configuration may negatively affect waste removal and nutrients supply. In agreement with these observations, in vitro studies demonstrated that prevailing low oxygen levels can impact on bovine granulosa cell, inducing their early luteinization with down-regulation of FSH signaling, cell proliferation and steroidogenesis beside up-regulation of angiogenesis, glucose metabolism, and inflammatory processes [79].
Finally, the alginate-filled configuration seems to have been the best analyzed condition, since the oocytes that matured in this condition showed an improvement in their oxidative/bioenergetic status, even at a maturation rate comparable with controls. The higher COC proximity, concave support curvature, 3D environment provided by the gel [80], and the absence of direct flow and basal oxygen levels may have recreated an environment similar to that of the ovarian follicle [80].
In conclusion, the obtained results demonstrate that, in the sheep model, mIVM with biocompatible supports protects the COCs from direct medium flow and elevated oxygen concentrations while keeping optimal nutrient supply and waste removal. This result is a prelude to an improvement of the developmental competence of oocytes matured under mIVM. This culture method could be applied to oocytes in suboptimal conditions, such as those recovered from ovine breeds under risk of genetic erosion, from other animal endangered species, from aged donors, or to oocytes matured after cryopreservation.

Conflicts of Interest:
The authors declare no conflict of interest.


Identification of Heme Oxygenase-1 as a Putative DNA-Binding Protein

Heme oxygenase-1 (HO-1) is a rate-limiting enzyme in degrading heme into biliverdin and iron. HO-1 can also enter the nucleus and regulate gene transcription independent of its enzymatic activity. Whether HO-1 can alter gene expression through direct binding to target DNA remains unclear. Here, we performed HO-1 CHIP-seq and then employed 3D structural modeling to reveal putative HO-1 DNA binding domains. We identified three probable DNA binding domains on HO-1. Using the Proteinarium, we identified several genes as the most highly connected nodes in the interactome among the HO-1 gene binding targets. We further demonstrated that HO-1 modulates the expression of these key genes using Hmox1 deficient cells. Finally, mutation of four conserved amino acids (E215, I211, E201, and Q27) within HO-1 DNA binding domain 1 significantly increased expression of Gtpbp3 and Eif1 genes that were identified within the top 10 binding hits normalized by gene length predicted to bind this domain. Based on these data, we conclude that HO-1 protein is a putative DNA binding protein, and regulates targeted gene expression. This provides the foundation for developing specific inhibitors or activators targeting HO-1 DNA binding domains to modulate targeted gene expression and corresponding cellular function.

Introduction
Heme oxygenase-1 (HO-1) is the rate-limiting enzyme of heme metabolism, which catalyzes the degradation of the heme into biliverdin, iron and carbon monoxide. It is a stress-inducible protein and one of the most ubiquitous and highly conserved antioxidant enzymes. A myriad of functions have been attributed to HO-1 as it plays important roles in maintaining cell homeostasis, tumorigenesis, immune function, reducing oxidative damage, and modulating both cell proliferation and apoptosis via regulation of intracellular levels of heme and heme metabolites [1][2][3][4]. Interestingly, HO-1 is also able to modulate cell function independently of its enzymatic activity by binding with itself, its isoenzyme HO-2, and caveolin-1 through protein-protein interactions [5]. As an example, HO-1 activity is negatively modulated by caveolin-1, by binding to a specific five-amino-acid long sequence in caveolin-1. Hyperoxia was shown to induce the physical interaction between caveolin-1 and HO-1 in fibroblasts, which reduced HO enzymatic activity [6]. As a compensatory response, oxidative stress led to activation of the antioxidant response element nuclear factor erythroid 2-related factor 2 (Nrf2), which induces HO-1, resulting humidified chamber. The cells for ChIP assay were grown in 15 cm 2 culture plates and collected at 75% confluency.

CRISPR/Cas9 Deficient of HO-1 in MLE-12 Cells
A two-step CRISPR system (Sigma-Aldrich, St. Louis, MO, USA) was used as per the manufacturer's instructions. The Cas9-blasticidin (LVCAS9BST) lentivirus was transduced into MLE-12 cells followed by blasticidin selection. Cells were then transduced with a noncoding-gRNA lentivirus and a HO-1-gRNA (TATGTAAAGCGTCTCCACG) puromycin + lentivirus. Cells were then selected with puromycin and single-sorted into 96-well plates using flow cytometry. Deletion efficiency was determined by measuring HO-1 protein levels using Western blot.

Transcriptional Analysis
qRT-PCR assays were performed using the ThermoFisher Taqman qPCR master mix and Taqman fluorescent probes. The qPCR assays were done with 18s (Taqman Hs99999) as a housekeeping control in an ABI5000 instrument. Probes used were Ubc (Mm01198158_m1), Dnajb1

ChIP-Seq Assay
The ChIP-Seq experiment was adapted from previously described [21]. MLE-12 cells were incubated with 1% formaldehyde for 10 min for crosslinking. After washing with wash buffer (10 mM HEPES, 10 mM EDTA, 0.5 mM EGTA, 0.25% Triton X-100, 0.2 mM PMSF, and 1 × protease inhibitors), the cells were counted and aliquoted into 100 × 10 6 cells per sample. Cells were then sonicated to shear chromatin using lysis buffer (1% SDS, 20 mM HEPES, 10 mM EDTA, and 50 mM NaCl) in a Covaris S220 High Performance Ultrasonicator. To confirm that the DNA shearing process was successful, the DNA was incubated with RNAseA (Thermo Scientific, EN0531, Waltham, MA, USA) and reversecrosslinked overnight at 65 • C followed by proteinase K incubation (Thermo Scientific, EO0491), phenol chloroform extraction, and ethanol precipitation. The clean DNA was run on a fragment analyzer to confirm desired shearing around 200 bp in length. The sheared chromatin was then incubated with HO-1 antibody (Abcam, ab13243, Waltham, MA, USA) or IgG (Abcam, ab171870) overnight with dynabeads A/G beads (Thermo Scientific, 78609). The beads were magnetically isolated and washed with RIPA buffer (50 mM TRIS, 1% NP-40, 2 mM EDTA, 0.1% sodium deoxycholate, 0.1% SDS, 300 mM NaCl) and DNA was isolated as per manufacturer's instructions. The ChIP-seq library was generated by Genewiz (Genewiz, South Plainsfield, NJ, USA) using Truseq ChIP library generation kit and sequenced with Illumina Hi-Seq 2 × 150 bp single index per lane with an average of 58 million reads per sample ( Figure 1A,B).
Further downstream analysis was performed in R 3.6.3 (R Core Team (2020). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria). Non-specific peaks from the IgG libraries were filtered out using findOverlapsOfPeaks from the ChIPPeakAnno (V 3.20.1) package [25,26]. The resulting peaks were annotated using annotatePeak from the ChIPseeker (V 1.22.1) package [27]. The distribution of peaks across promoters, exons, introns, UTRs, and distal or intergenic regions were visualized using ggplot2 (V 3.3.0 H. Wickham. ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag New York, 2016.) and the annoStat slot from the anno-tatePeak results. The distribution of peaks relative to transcription start sites were visualized using ChIPseeker plotDistToTSS function [27].
Further downstream analysis was performed in R 3.6.3 (R Core Team (2020). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria). Non-specific peaks from the IgG libraries were filtered out using findOverlapsOfPeaks from the ChIPPeakAnno (V 3.20.1) package [25,26]. The resulting peaks were annotated using annotatePeak from the ChIPseeker (V 1.22.1) package [27]. The distribution of peaks across promoters, exons, introns, UTRs, and distal or intergenic regions were visualized using ggplot2 (V 3.3.0 H. Wickham. ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag New York, 2016.) and the annoStat slot from the annotate-Peak results. The distribution of peaks relative to transcription start sites were visualized using ChIPseeker plotDistToTSS function [27].
All of the code for computational analysis of ChIP data (alignment, read and alignment QC, read trimming, and downstream analyses in R) is available on GitHub (https://github. com/compbiocore/Scaffa, accessed on 1 September 2022).
Hit mapped reads were further analyzed for their location in the genome using inhouse scripts and using data from the UCSC human genome browser [28], gene ontology was observed using gclusterProfiler to run enrichment analysis (enrichKEGG function), and common motifs of the reads were found de novo using the MEME suite (MEME Suite Software) [29]. Statistical analyses such as t-tests and one-way ANOVA were ran using Prism 9.4.0 (GraphPad Software, San Diego, CA, USA).

Performing Protein-Protein Interaction (PPI) Networks (Proteinarium)
We ran Proteinarium [18] on the genes shown to bind HO-1 by ChIP-Seq and on the genes identified through our BLAST search [30] using the results of the MEME motif search [31]. Proteinarium is a PPI network visualization and analysis tool. Proteinarium exclusively accepts human gene identifiers so we converted mouse gene names to HGNC (HUGO Gene nomenclature committee) symbols prior to analysis [32]. Out of 118 genes identified by ChIP-Seq, we excluded 48 genes which did not carry either official names or interaction information. Proteinarium was run using default parameters and by setting the maximum path length to two nodes of separation. Here, the maximum path length parameter defines the number of imputed proteins which may link seed genes using experimentally validated protein-protein interaction data.

Functional Enrichment
We performed Gene Ontology (GO) Term [33] and Kyoto Encyclopedia of Genes and Genomes (KEGG) [34] pathway enrichments using the Functional Annotation Tool provided by the Database for Annotation, Visualization and Integrated Discovery (DAVID) v6.8 [35,36]. Bubble charts of significantly enriched Biological Process GO Terms and their associated False Discovery Rate corrected p values were plotted using the ggplot2 v3.3.2 R package.

Predictive Modeling of Apo HO-1 dsDNA Binding Complex Structures
We developed a novel workflow for predictive structural modeling of the dsDNA binding complex of human apo HO-1 and characterization of genomic binding targets of human apo HO-1 using ModelX [19], FoldX [37], UCSF Chimera [38], Linux, and R scripting. The complete workflow consists of three main steps: (1) all possible dsDNA binding sites and dsDNA backbone conformations on HO-1 were identified using ModelX software; (2) the predicted set of all possible dsDNA backbone binding fragments to a set of continuous and highly energetically favorable dsDNA-peptide complexes in three distinct locations on the 3D structure of HO-1 were refined using FoldX software; and (3) the specific dsDNA nucleotide sequences within each dsDNA backbones which were most energetically favorable to bind the peptide residues in each putative DNA binding domain were identified using FoldX software. Finally, we showed the molecular binding interactions between each of the three DNA binding domains and their compatible double stranded DNA molecules using the visualization program Chimera. We have described our workflow in detail below.

ModelX-dsDNA-Protein Docking
We used the ModelX tool suite for biomolecular modeling to identify short dsDNA backbone structures which are compatible binding partners with residues in a query Protein Data Bank (PDB) protein file. We selected the human heme-free apo HO-1 protein PDB structure which we extracted from the PBD file (PDB ID = 1N16) as a query protein structure since it has been previously shown that the heme binding C-terminus domain is cleaved from HO-1 before it's transportation to the nucleus [12]. We used selective filters to find many short dsDNA helix backbone fragments which are predicted to bind HO-1 with high confidence and assembled a cloud of all predicted fragments positioned around the HO-1 structure in PDB format. We then searched for continuous stretches of short dsDNA backbones in the model to identify extended fragments in order to further refine our set of predicted dsDNA binding backbones. By identifying backbone sequences with greater lengths, we were able to identify binding sequences contained within the modeled dsDNA helix backbones with higher specificity since they provide more points of comparison.

Identification of Compatible dsDNA Binding Targets of Apo HO-1
We performed a series of energy calculations using the DNAScan tool from the FoldX Suite [21] to predict the most probable nucleotide sequence contained within each extended dsDNA backbone fragment. Using the energy file output of DNAScan, we assembled energy matrices containing the change in Gibbs Free Energy (∆∆G) of the HO-1 structure incurred by mutating each base along each of the predicted binding fragments to each of the four nucleotides. We used these energy matrices to identify the most probable dsDNA sequences contained in each of the dsDNA backbone fragments by selecting the nucleotides in each fragment position which produced the most stable binding complex structure with the lowest free energy.
We used the DNAContact FoldX tool [21] to identify the positions along each dsDNA fragment which make direct contact with one or more residues on the HO-1 structure. We only considered fragments which had at least three DNA contact positions for further analysis.
To identify compatible dsDNA genomic binding targets of HO-1, we used these newly identified energetically favorable binding sequences to thread the genomic sequence of our genes of interest. We only considered positions within the dsDNA fragments which were shown to directly bind HO-1 by the results from FoldX DNAContact. Since our calculated energetically favorable binding sequences are double stranded, we simultaneously threaded both the genomic sequence and its complementary sequence to identify compatible binding sequences.
We identified all compatible genomic binding positions in our target genes and recorded the coordinates as well as total number of binding sites per gene. We identified genes which are predicted to have the highest binding affinity to HO-1 by normalizing the total number of compatible HO-1 binding sites on each gene by the gene length. We characterized the functions of the genes with the greatest number of compatible binding sites using DAVID functional enrichment tools. We performed sequence threading using the sequences of the genes in the gene lists produced by CHIP-seq and by our subsequent motif search and have summarized this data in several tables (Supplementary Tables S1 and S7).

Visualization of HO-1 Protein and dsDNA Binding Fragment PDB Data
We used UCSF Chimera program to perform all visualizations of dsDNA fragments and human apo HO-1 structural data in PDB format [21].

Transfection of Mutated HO-1 Plasmids into Hmox1 Deficient Cells
Four conserved amino acids, including E215, I211, E201, and Q27, in HO-1 DBD1 were mutated to alanine, which was validated by sequencing and cloned in the pcDNA3.1 vector. The wild-type or mutated HO-1 plasmids (2 µg) were transfected into Hmox1 deficient MLE-12 cells for 24 h by using lipofectamine 3000 according to the manufacturer's protocol. The cell pellets were collected for measuring expression of Fzd2, Gtpbp3, Rpl36, and Eif1 genes using the qRT-PCR.

Statistical Analysis
Experiments were repeated biologically three times, with three technical replicates for quantitative qRT-PCR assay. Data are expressed as mean ± standard error of the mean (SEM). The unpaired t-test was used for detecting statistical significance of the differences between means of two groups if the data are normally distributed using Graph Pad Prism 9 software. Welch-corrected t test was used if the data are not normally distributed. For multiple comparisons, the statistical significance of the differences was evaluated by using one-way ANOVA followed by Tukey's post-test to specifically compare indicated groups. Statistical significance was considered when the p value was less than 0.05.

CHIP-Seq Reveals That HO-1 Binds to DNA
ChIP-Seq is a robust method used to analyze protein interactions with DNA. Thus, we performed a ChIP-Seq in mouse lung epithelial MLE-12 cells to determine whether HO-1 could be a DNA-binding protein. The experiment yielded over 50 million 250 pairedend reads per sample. The data went through quality control, alignment, and further downstream analysis (Supplementary Tables S1 and S2, Supplementary Figures S1-S3). Compared to IgG control, hundreds of hits on 118 separate genes were observed in the HO-1 ChIP assay (Supplementary Table S7). These aligned to less than 1 kb from the promoter (33.97% of hits) and within the distal intergenic regions (37.8% of hits) ( Figure 1A-C). The KEGG gene ontology analysis of the DNA hits revealed that the genes are part of pathways related to cellular responses to external stimuli, cellular responses to stress, cellular senescence, DNA replication pre-initiation, oxidative stress-induced senescence, programmed cell death, G1/S transition, assembly of the pre-replicative complex, metalloproteases DUBs, and apoptosis-induced DNA fragmentation (Supplementary Figure S4). Table 1 showed the top 10% of the 118 genes with highest statistical significance as well as a known relationship of the genes to HO-1 in the literature. The sequences from all 118 gene targets were submitted to motif analysis yielding 10 motifs that are common across subsets of our hits (Figure 2A). Although none of these motifs were a perfect match to known motif libraries, variations of these motifs were found to align to 130 genes including 45 genes that had multiple motifs within each gene. From these motifs we also ran function analysis leading to the molecular functions, biological processes, and cellular components in which the genes are known to play a role ( Figure 2B).
These results suggest that HO-1 protein putatively binds to DNA with a preference for specific motifs sequences.
The sequences from all 118 gene targets were submitted to motif analysis yielding 10 motifs that are common across subsets of our hits (Figure 2A). Although none of these motifs were a perfect match to known motif libraries, variations of these motifs were found to align to 130 genes including 45 genes that had multiple motifs within each gene. From these motifs we also ran function analysis leading to the molecular functions, biological processes, and cellular components in which the genes are known to play a role ( Figure 2B). These results suggest that HO-1 protein putatively binds to DNA with a preference for specific motifs sequences.

PPI Network Analysis Reveals Multiple Interactions
To get additional confirmation of binding interactions and binding sites, we analyzed the PPI network of the gene set produced from the ChIP-Seq results and produced a network visualization of the interactions by using the Proteinarium. We ran Proteinarium with the list of 70 human gene identifiers (ids) and generated a protein-protein interaction network that contains 34 genes. Of the 70 genes in the input list, 15 were found to interact in a network with 19 imputed genes based upon a maximum imputed path length of 2 nodes of separation. Using Proteinarium, we identified Uba52, Ubc, Hist2h2ab, Hist2h2ac, Ddx5, Eef1a1, Rpl3, Rpl9, Rpl27, Rpl36, and Hnrnpa3 as the most highly connected nodes in the interactome among the HO-1 gene binding targets identified by CHIP-seq ( Figure  3A). These results suggest that there are multiple genes that interact with HO-1 protein.

PPI Network Analysis Reveals Multiple Interactions
To get additional confirmation of binding interactions and binding sites, we analyzed the PPI network of the gene set produced from the ChIP-Seq results and produced a network visualization of the interactions by using the Proteinarium. We ran Proteinarium with the list of 70 human gene identifiers (ids) and generated a protein-protein interaction network that contains 34 genes. Of the 70 genes in the input list, 15 were found to interact in a network with 19 imputed genes based upon a maximum imputed path length of 2 nodes of separation. Using Proteinarium, we identified Uba52, Ubc, Hist2h2ab, Hist2h2ac, Ddx5, Eef1a1, Rpl3, Rpl9, Rpl27, Rpl36, and Hnrnpa3 as the most highly connected nodes in the interactome among the HO-1 gene binding targets identified by CHIP-seq ( Figure 3A). These results suggest that there are multiple genes that interact with HO-1 protein.

Functional Enrichment of Protein-Protein Interaction Network Reveals Many Biological Processes Are Impacted by HO-1
Functional enrichment analysis of the 34 genes found to exist in the protein-protein interaction network yielded significantly enriched GO and KEGG pathway terms using a false discovery rate (FDR) corrected p-value threshold of 0.05 or lower. Top significantly enriched biological process GO terms among genes in the protein-protein interaction network included mRNA splicing via the spliceosome, gene expression, nuclear-transcribed mRNA catabolic processes, SRP-dependent cotranslational protein targeting, rRNA processing, DNA damage response, signal transduction by p53 class mediator resulting in cell cycle arrest, translational initiation, and fibroblast growth receptor signaling ( Figure 3B). Top significantly enriched cellular component GO terms included the nucleoplasm, intracellular ribonucleoprotein complex, catalytic step 2 spliceosome, nucleus, nucleolus, nucleosome, cytosolic large ribosomal subunit, extracellular exosome, and nuclear chromatic. Two significantly enriched molecular function GO Terms were identified to be included poly(A) RNA binding and RNA binding. Additionally, the spliceosome KEGG pathway was identified to be significantly enriched. Altogether, these results indicate that HO-1 modulates gene expression involved in many biological processes.

Functional Enrichment of Protein-Protein Interaction Network Reveals Many Biological Processes Are Impacted by HO-1
Functional enrichment analysis of the 34 genes found to exist in the protein-protein interaction network yielded significantly enriched GO and KEGG pathway terms using a false discovery rate (FDR) corrected p-value threshold of 0.05 or lower. Top significantly enriched biological process GO terms among genes in the protein-protein interaction network included mRNA splicing via the spliceosome, gene expression, nuclear-transcribed mRNA catabolic processes, SRP-dependent cotranslational protein targeting, rRNA processing, DNA damage response, signal transduction by p53 class mediator resulting in cell cycle arrest, translational initiation, and fibroblast growth receptor signaling ( Figure  3B). Top significantly enriched cellular component GO terms included the nucleoplasm, intracellular ribonucleoprotein complex, catalytic step 2 spliceosome, nucleus, nucleolus, nucleosome, cytosolic large ribosomal subunit, extracellular exosome, and nuclear chromatic. Two significantly enriched molecular function GO Terms were identified to be included poly(A) RNA binding and RNA binding. Additionally, the spliceosome KEGG pathway was identified to be significantly enriched. Altogether, these results indicate that HO-1 modulates gene expression involved in many biological processes.

Predictive Modeling of Apo HO-1 dsDNA Binding Complex Structures Reveals Several HO-1 DNA Binding Domains
To identify the predicted dsDNA binding domains on human apo heme-free HO-1, we used the workflow powered by ModelX and FoldX software. We identified four predicted dsDNA binding domains on human apo heme-free HO-1. We excluded one binding domain for which the predicted dsDNA fragment had only one DNA contact position. The other three binding sites were included for consideration as probable DNA binding domains on HO-1. Here, we have named these predicted binding domains HO-1 DBD-1, DBD-2, and DBD-3 ( Figure 4). We presented our predictive characterization of these three human heme-free apo HO-1 DNA-binding domains ( Figure 5). Summary and detailed information on the number of binding protein residues and nucleotides, specific compatible nucleotides and amino acids for each of the three DNA binding domains are described in detail in Tables 2 and 3. Three probable DNA binding domains of HO-1 are identified. domains on HO-1. Here, we have named these predicted binding domains HO-1 DBD-1, DBD-2, and DBD-3 ( Figure 4). We presented our predictive characterization of these three human heme-free apo HO-1 DNA-binding domains ( Figure 5). Summary and detailed information on the number of binding protein residues and nucleotides, specific compatible nucleotides and amino acids for each of the three DNA binding domains are described in detail in Tables 2 and 3. Three probable DNA binding domains of HO-1 are identified.

Identification of Compatible dsDNA Binding Targets of Apo HO-1 among CHIP-Seq Genes
As shown in the link (https://docs.google.com/spreadsheets/d/1ewOX-wBOIg69V1 DAs8Ua5HCDXALqsPFr/edit?usp=sharing&ouid=117066342160448117876&rtpof=true& sd=true, accessed on 1 September 2022), we identified targeted CHIP-seq genes that bind to HO-1 DNA-binding domains. The murine targets were converted to human targets since there is no published murine apo HO-1 protein structure. We describe the top 10 genes predicted to bind each HO-1 DNA-binding domain. Several of the genes with the highest HO-1 binding affinity, including UBA52, Rpl27, and Rpl36, are among the most connected nodes in the PPI network generated with Proteinarium.

HO-1 DNA Binding Domains DBD-1 DBD-2 DBD-3
Predicted binding dsDNA fragment length (number of nucleotides) 10 10 10 Protein residue contacts 9 7 8 DNA contacts 7 6 5 Here, we report the fragment length (in number of nucleotides) of each predicted dsDNA binding fragment predicted by ModelX, the number of HO-1 protein residues contacting the dsDNA fragments, and the number of DNA fragment positions which contact with protein residues for each of the newly described three HO-1 DNA binding domains.   Here, we report the fragment length (in number of nucleotides) of each predicted dsDNA binding fragment predicted by ModelX, the number of HO-1 protein residues contacting the dsDNA fragments, and the number of DNA fragment positions which contact with protein residues for each of the newly described three HO-1 DNA binding domains.

HO-1 DNA-Binding Domain 1 Targets among CHIP-Seq Genes
We found that 62 of 71 human genes which were found to bind HO-1 by ChIP-Seq and which are conserved between human and mice contain at least 1 dsDNA binding sequence which is compatible with the HO-1 DBD-1 binding domain. Of these 62 genes, Hmgb2, Fzd2, Gtpbp3, Rpl36, Ankle1, Mrpl12, Npm3, Particl, Eif1, and Ranbp1 were the top 10 genes with the greatest proportion of HO-1 DBD-1 compatible binding sites normalized by gene length.

HO-1 DNA-Binding Domain 2 Targets among CHIP-Seq Genes
We found that 65 of 71 human genes which were found to bind HO-1 by ChIP-Seq and which are conserved between human and mice contain at least 1 dsDNA binding sequence which is compatible with the HO-1 DBD-2 binding domain. Of these 65 genes, Rpl27, Hmgb2, Ankle1, Mat2a, Dnajb1, Mrpl12, Npm3, Nfe2l3, Orc6, and Ranbp1 were the top 10 genes with the greatest proportion of HO-1 DBD-2 compatible binding sites normalized by gene length.

HO-1 DNA-Binding Domain 3 Targets among CHIP-Seq gENE
We found that 68 of 71 human genes which were found to bind HO-1 by ChIP-Seq and which are conserved between human and mice contain at least 1 dsDNA binding sequence which is compatible with the HO-1DBD-3 binding domain. Of these 68 genes, Rn7sk, Npm3, Particl, Malat1, Etfrf1, Hmgb2, Sart3, Nkain2, Fut10, and Uba52 were the top 10 genes with the greatest proportion of HO-1 DBD-3 compatible binding sites normalized by gene length.

Identification of Compatible dsDNA Binding Targets of Apo HO-1 among Motif BLAST Hit Genes
As shown in the link (https://docs.google.com/spreadsheets/d/1ux5PoBpGpWJ5 iaeuCpUf-PkFrrGGkE7x/edit?usp=sharing&ouid=117066342160448117876&rtpof=true& sd=true, accessed on 1 September 2022), we identified targeted Motif BLAST Hit genes that bind to HO-1 DNA-binding domains. The murine targets were also converted to human targets since there is no published murine apo HO-1 protein structure.

HO-1 DNA-Binding Domain 1 Targets among Motif BLAST Hit Genes
We found that 99 of 100 human genes which were found to have close homology with DNA binding motifs identified from genes found to bind HO-1 by ChIP-Seq contain at least 1 dsDNA binding sequence which is compatible with the HO-1 DBD-1 binding domain. Of these 99 genes, Foxq1, Etv2, Ascl2, Foxc1, Foxc2, Sp5, Barx1, Sox4, and Maz were the top 10 genes with the greatest proportion of HO-1 DBD-1 compatible binding sites normalized by gene length.

HO-1 DNA-Binding Domain 2 Targets among Motif BLAST Hit Genes
We found that 98 of 100 human genes which were found to have close homology with DNA binding motifs identified from genes found to bind HO-1 by ChIP-Seq contain at least 1 dsDNA binding sequence which is compatible with the HO-1 DBD-2 binding domain. Of these 99 genes, Sp5, Maz, Barx1, Foxa2, Foxa1, Ascl2, Lyl1, Rxra, Srf, and Evx1 were the top 10 genes with the greatest proportion of HO-1DBD-2 compatible binding sites normalized by gene length.

HO-1 DNA-Binding Domain 3 Targets among Motif BLAST Hit Genes
We found that 99 of 100 human genes, which were found to have close homology with DNA binding motifs identified from genes found to bind HO-1 by ChIP-Seq, contain at least 1 dsDNA binding sequence compatible with the HO-1 DBD-3 binding domain. Of these 99 genes, Sry, Nr1d1, Tbx21, E2f4, Gli1, Sox2, Foxq1, Tbx3, Nanog, and Foxj2 were the top 10 genes with the greatest proportion of HO-1 DBD-3 compatible binding sites normalized by gene length.

Functional Enrichment of HO-1 Binding Gene Targets Reveals a Role in Regulation of Transcription Factor Regulation, Cell Proliferation, and Tissue Development
To understand the functional impact of predicted HO-1 binding on regulatory domains of the genes found to bind HO-1 by ChIP-Seq, we performed a functional enrichment analysis on the top 50% of genes which had a higher proportion of HO-1 binding sequences after normalization by gene length. We found that the top 50% of gene targets for all three DBDs had commonly shared significantly enriched biological functions. The top significantly enriched biological functions shared by the gene targets of all three HO-1 DBDs involved DNA transcription, regulation of transcription, and tissue development including hemopoiesis, regulation of cell proliferation and organogenesis (Supplementary Figure S4).

Hmox1 Deficient MLE-12 Cells Show Dysregulation of Multiple Target Genes
In order to determine whether HO-1 had a direct effect on the transcription of identified genes, we generated a Hmox1 deficient line using the CRISPR-Cas9 system. Using qRT-PCR, we showed that expression of several genes identified as binding to HO-1, namely DnaJ heat shock protein family Hsp40/Dnajb1 Member B1 (Hsp40), eukaryotic translation elongation factor 1 alpha 1 (Eefa1), ribosomal protein L27 (Rpl27), ubiquitin C (Ubc), and high-mobility group protein B2 (Hmgb2), was significantly downregulated in Hmox1 deficient MLE-12 cells (Figure 6, Supplementary Figure S5). In contrast, the expression of nuclear receptor subfamily 1 group d member 1 (Nr1d1), metastasis associated lung adenocarcinoma transcript 1 (Malat1), and forkhead box a1 (Foxa1) genes was significantly upregulated in Hmox1 deficient cells ( Figure 6). There were no changes in E2f transcription factor 1 (E2f1), eukaryotic translation initiation factor 1 (Eif1), heterogeneous nuclear ribonucleoprotein a2/b11 (Hnrnpa2b), or ubiquitin a-52 residue ribosomal protein fusion product 1 (Uba52) genes between Hmox1 deficient and WT cells ( Figure 6). This suggests that HO-1 modulates gene expression in a gene-specific manner, and not all genes that are predicted to bind to HO-1 protein are transcriptionally regulated under these conditions. suggests that HO-1 modulates gene expression in a gene-specific manner, and not all genes that are predicted to bind to HO-1 protein are transcriptionally regulated under these conditions.

Mutation of Amino Acids in HO-1 DBD1 Alters Target Gene Expression
To further validate HO-1 as a DNA binding protein, we constructed plasmids of wild type and mutated HO-1 ( Figure 7A,B). In mutated plasmids, four amino acids (E215, I211, E201, and Q27) within the HO-1 DBD1 were mutated to alanine. These were chosen since they are conserved between mouse and human. We then transfected these plasmids into Hmox1 deficient MLE-12 cells and measured expression of Fzd2, Gtpbp3, Rpl36, and Eif1 genes. We chose these genes, as they were identified within the top 10 binding hits normalized by gene length predicted to bind DBD-1, but are not in the top 10 binding hits for DBD-2 or DBD-3. As shown in Figure 7C, mutation of these amino acids within the HO-1 DBD-1 did not affect Hmox1 gene expression, but significantly increased the expression of Gtpbp3 and Eif1 genes. These results suggest that amino acids E215, I211, E201, and Q27 within the DBD1 play important roles in binding and regulating expression of target genes.

Mutation of Amino Acids in HO-1 DBD1 Alters Target Gene Expression
To further validate HO-1 as a DNA binding protein, we constructed plasmids of wild type and mutated HO-1 ( Figure 7A,B). In mutated plasmids, four amino acids (E215, I211, E201, and Q27) within the HO-1 DBD1 were mutated to alanine. These were chosen since they are conserved between mouse and human. We then transfected these plasmids into Hmox1 deficient MLE-12 cells and measured expression of Fzd2, Gtpbp3, Rpl36, and Eif1 genes. We chose these genes, as they were identified within the top 10 binding hits normalized by gene length predicted to bind DBD-1, but are not in the top 10 binding hits for DBD-2 or DBD-3. As shown in Figure 7C, mutation of these amino acids within the HO-1 DBD-1 did not affect Hmox1 gene expression, but significantly increased the expression of Gtpbp3 and Eif1 genes. These results suggest that amino acids E215, I211, E201, and Q27 within the DBD1 play important roles in binding and regulating expression of target genes.

Discussion
Previous studies have shown that HO-1 protein can enter the nucleus and bind to nuclear proteins, which regulate transcriptional factors including AP1, AP2, Brn3, Stat1, Stat3, and Stat4 [11,12,14,39]. In addition, the DNA binding activity of NF-κB and stimulatory protein (SP)-1 are repressed when both nuclear and cytoplasmic HO-1 protein is overexpressed 21. In addition, during oxidative stress, nuclear HO-1 binds to Nrf2 to prevent its degradation, thereby prolonging its action [40]. It is unclear whether HO-1 can directly bind to DNA and modulate gene expression. In the present study, using Chipseq, we demonstrated, for the first time, that HO-1 is a putative DNA binding protein and identified three probable HO-1 DNA binding domains using 3D structural modeling. By using protein-protein interaction networks, we identified Uba52, Ubc, Hist2h2ab, Hist2h2ac, Ddx5, Eef1a1, Rpl3, Rpl9, Rpl27, Rpl36, and Hnrnpa3 as the most highly connected nodes in the interactome among the HO-1 gene binding targets identified by CHIPseq. We validated the impact of HO-1 mutations within the DBD-1 and disruption on the expression of the target genes. Therefore, HO-1 protein putatively binds to DNA, which may regulate target gene expression. The expression of HO-1 is low under normal physiological conditions. Further studies are warranted to investigate the effect of HO-1 on target genes under oxidative stress.
Unlike a DNA electrophoretic mobility shift assay, ChIP-seq allows for the characterization of a protein that binds to DNA as well as the sequences of the oligos from the immunoprecipitated protein-DNA complex [41]. Using this technology, we showed that HO-1 bound to specific DNA sequences. Future studies using surface plasmon resonance and isothermal titration calorimetry as complementary methods would further validate the interaction of HO-1 with DNA. One third of the DNA hits were within promoter regions, suggesting that this could alter transcriptional regulation. In addition, another third of the binding sites found within distal intergenic regions could still have an effect on gene transcription and processing. Transcription factors bind to enhancer or promoter regions

Discussion
Previous studies have shown that HO-1 protein can enter the nucleus and bind to nuclear proteins, which regulate transcriptional factors including AP1, AP2, Brn3, Stat1, Stat3, and Stat4 [11,12,14,39]. In addition, the DNA binding activity of NF-κB and stimulatory protein (SP)-1 are repressed when both nuclear and cytoplasmic HO-1 protein is overexpressed 21. In addition, during oxidative stress, nuclear HO-1 binds to Nrf2 to prevent its degradation, thereby prolonging its action [40]. It is unclear whether HO-1 can directly bind to DNA and modulate gene expression. In the present study, using Chip-seq, we demonstrated, for the first time, that HO-1 is a putative DNA binding protein and identified three probable HO-1 DNA binding domains using 3D structural modeling. By using protein-protein interaction networks, we identified Uba52, Ubc, Hist2h2ab, Hist2h2ac, Ddx5, Eef1a1, Rpl3, Rpl9, Rpl27, Rpl36, and Hnrnpa3 as the most highly connected nodes in the interactome among the HO-1 gene binding targets identified by CHIP-seq. We validated the impact of HO-1 mutations within the DBD-1 and disruption on the expression of the target genes. Therefore, HO-1 protein putatively binds to DNA, which may regulate target gene expression. The expression of HO-1 is low under normal physiological conditions. Further studies are warranted to investigate the effect of HO-1 on target genes under oxidative stress.
Unlike a DNA electrophoretic mobility shift assay, ChIP-seq allows for the characterization of a protein that binds to DNA as well as the sequences of the oligos from the immunoprecipitated protein-DNA complex [41]. Using this technology, we showed that HO-1 bound to specific DNA sequences. Future studies using surface plasmon resonance and isothermal titration calorimetry as complementary methods would further validate the interaction of HO-1 with DNA. One third of the DNA hits were within promoter regions, suggesting that this could alter transcriptional regulation. In addition, another third of the binding sites found within distal intergenic regions could still have an effect on gene transcription and processing. Transcription factors bind to enhancer or promoter regions to cause an effect on target gene transcription. It is unclear whether HO-1 interacts directly with transcription factors, the spliceosome, coactivators or corepressors, which co-modulates targeted gene expression. It will also be useful to verify our findings in other cell types to determine whether the binding of HO-1 protein to DNA is ubiquitous and not cell-specific. A luciferase reporter assay is warranted to establish a functional connection between the presence of HO-1 protein and the amount of gene product synthesized.
The murine apo HO-1 protein structure is not available in the Protein Data Bank. Many domains of HO-1 protein, including the alpha helical fold, are highly conserved among bacteria, mice, and humans [16,42,43]. Thus, we used the human apo HO-1 structure for modeling with ModelX and FoldX software. The murine structure would have been preferred since the ChIP-Seq was done in murine cells. We identified three putative novel DNA binding domains on human apo HO-1 protein with nucleotide resolution of specific compatible double-stranded DNA (dsDNA) binding sequences. Crystallography could conclusively demonstrate the co-crystallization of the DNA targets found here.
Using the refined nucleotide sequences for each of the unique dsDNA fragments predicted to bind the three predicted DBDs, we were able to search the genes found to bind HO-1 by ChIP-Seq for compatible binding sites within each gene. We showed that most genes found to bind HO-1 by ChIP-Seq had at least one predicted HO-1 compatible binding sequence and identified genes which had the highest proportion of HO-1 compatible binding sequences normalized by their gene length. The most energetically favorable binding sequences for all three sites had a high proportion of GC content, which suggests that HO-1 DNA binding favors GC rich regions such as gene regulatory domains. These GC-rich sequences are often located in gene promoters and play a role in transcription initiation. For instance, 76% of human promoters are GC-rich and contain multiple binding sites of the transcriptional activator SP1 [44]. Whether these GC-rich regions are located in the promoters of HO-1 targeted genes remains unknown.
We observed that HO-1 DBD-2 overlaps with the heme-binding functional ring of HO-1 on peptide residue His25. It is possible that HO-1 His25 is not available for binding heme when dsDNA is bound to DBD-2. This could explain the loss of HO-1 enzymatic activity when HO-1 enters the nucleus. It remains unclear whether heme is competing for the binding of HO-1 to DNA, although it is possible since heme has been found in small quantities within the nucleus [45]. DBD-1 contains Gln27 on chain B, but the heme binding domain is on chain A (opposite side of the protein). However, DBD-2 is on chain A and its peptides surround the heme binding functional ring on His25. Mutation of His25 significantly reduces HO-1 activity [46]. It remains unclear whether DBD-1 and DBD-2 domains are essential for HO-1 protein folding or enzymatic activity. This could be validated using our HO-1 DBD-1 mutated cells.
HO-1 can undergo post-translational modifications (PTMs), including phosphorylation and acetylation, which affects its localization and function [14,[47][48][49]. We observed that DNA binding domains contain or are very proximal to amino acids known to have PTMs in HO-1. This may affect or allow binding of HO-1 protein to targeted genes. For instance, DBD-1 is close to Tyr114 which can be phosphorylated and to Arg196 which can be acetylated. DBD-1 also contains the binding site of a known HO-1 inhibitor QC-308 [50]. Whether QC-308 causes inhibition of HO-1 by preventing HO-1 protein-DNA binding remains to be explored. While DBD-2 is likely to be regulated by ubiquitination at Lys22, other Lys residues including Lys18 and Lys 20 may promote ubiquitination or acetylation of HO-1 protein. Nevertheless, many residues in the DBDs have no PTMs [51]. Further study is required to determine whether mutation of these and other residues affects the binding of HO-1 protein to DNA. Importantly, screening drugs that only affect the DNA-binding function of HO-1 is an exciting future direction and may have important implications for health.
Of the gene targets explored, HO-1 deficiency caused transcriptional changes on Dnajb1, Eef1a1, Rpl27, Ubc, Hmgb2, Malat1, Foxa1, and Nr1d1. This suggests that HO-1 may at least partially regulate these genes through DNA transcription. These genes are involved in DNA flexibility, translation initiation, protein folding, differentiation, proliferation, and circadian oscillation [37,[52][53][54][55], corroborating the multiple functions of HO-1. Nr1d1 is also a target gene of Nrf2 [56]. Further studies are required to determine whether HO-1 interacts with other transcription factors including Nrf2 which regulates expression of these genes. A previous study comparing Hmox1 deficient and WT induced pluripotent stem cells reported differential gene expression of over 15 enriched biological process GO Terms [57]. Our study showed congruence of six biological process terms. This includes proliferation, differentiation, embryogenesis, DNA repair and tumorigenesis. These results can be used as a guide for further investigation of these potential targets.

Conclusions
HO-1 protein preferentially binds to DNA at specific DNA binding sequences or motifs within three novel putative binding sites, including one within its enzymatic pocket. This strongly suggests that the function of HO-1 is not only enzymatic, and regulatory through protein-protein binding, but also regulatory through DNA binding and transcriptional regulation. Thus, our study provides the first characterization of HO-1 as a putative DNA-binding protein and begins the foundation for developing specific inhibitors or activators targeting HO-1 DNA binding domains to modulate targeted gene expression and corresponding cellular function.


s2fieldsofstudy: Biology Medicine Chemistry

text:
Antiviral effects of deoxynojirimycin (DNJ)-based iminosugars in dengue virus-infected primary dendritic cells

Dendritic cells (DCs) are important targets for dengue virus (DENV) infection and play a significant role in the early immune response. Antiviral effects of iminosugars against DENV in primary cells have been demonstrated previously in monocyte-derived macrophages (MDMΦs). Given the important role played by DCs in innate immune defense against DENV, the antiviral effects of three deoxynojirimycin (DNJ) derivatives (NN-DNJ, EOO-DNJ and 2THO-DNJ) and a deoxygalactonojirimycin (DGJ) negative control were evaluated in DENV-infected primary human monocyte-derived immature DCs (imDCs). DNJ- but not DGJ-derivatives elicited antiviral activity in DENV-infected imDCs, similar to that observed in MDMΦs. The DNJ-derivatives inhibited DENV secretion in a dose-dependent manner. Endoplasmic reticulum (ER) α-glucosidase I inhibition by DNJ-derived iminosugars, at concentrations of 3.16 μM, correlated with a reduction in the specific infectivity of virions that were still secreted, as well as a reduction in DENV-induced tumour necrosis factor alpha secretion. This suggests iminosugar-mediated ER α-glucosidase I inhibition may give rise to further benefits during DENV infection, beyond the reduction in viral secretion associated with ER α-glucosidase II inhibition.

Introduction
Dengue virus (DENV) infections cause significant morbidity and mortality with increasing global incidence (Bhatt et al., 2013). While there has been recent progress with vaccine development (Wilder-smith, 2020), the quest for an effective antiviral continues. Certain iminosugars have long been under investigation as antiviral agents against several enveloped viruses (Zitzmann et al., 1999;Wu et al., 2002;Dwek et al., 2002), and previous studies have demonstrated that deoxynojirimycin (DNJ)-derived iminosugars reduce DENV secretion in a variety of infected cultured cells (Wu et al., 2002;Warfield et al., 2015). We have shown that DNJ-derived iminosugars have anti-DENV effects in primary monocyte derived macrophages (MDMΦs) and inhibit pathogen-induced inflammatory responses (Sayce et al., 2021). Promising in vitro and in vivo data for iminosugar antiviral activity against DENV have paved the way for iminosugars to enter clinical trials.
DNJ-derived iminosugars are glucose mimetics which possess a nitrogen atom in place of the oxygen in the carbohydrate ring; these inhibit host α-glucosidases pivotal to endoplasmic reticulum (ER) glycoprotein folding quality control (Alonzi et al., 2017). The outer α-1, 2-linked glucose is trimmed by ER α-glucosidase I (αGluI) and the inner α-1,3-linked glucoses are removed by ER α-glucosidase II (αGluII) during N-linked glycan processing (Chang et al., 2013). Impaired glucose trimming due to iminosugar treatment results in misfolded proteins which are either degraded by ER-associated degradation or secreted with altered properties. Iminosugar treatment of virally infected cells thus disrupts proper virion morphogenesis (Alonzi et al., 2017). Though iminosugars can inhibit both α-glucosidases and glycolipid processing pathways, previous work has shown that inhibition of α-glucosidases underlies the antiviral effect of iminosugars in DENV infection in MDMΦs (Sayce et al., 2016;Miller et al., 2018). While iminosugars inhibit αGluII at lower concentrations and αGluI at higher concentrations of compound, inhibition of αGluII alone leads to measurable anti-flaviviral activity (Kiappes et al., 2018). However, the potentially enhanced effects of combined αGluI and αGluII inhibition seen with higher iminosugar concentrations have not been fully explored in DENV infection.
Dendritic cells (DCs) play a crucial role in combating infection by diverse pathogens (Cheong et al., 2011;Nakano et al., 2009;Schmid and Harris, 2014) and DCs are the first cells to encounter the virus following a mosquito bite in DENV infection (Marovich et al., 2001). Tumour necrosis factor alpha (TNFα), secreted by immune cells including DCs (Ho et al., 2001), is an important mediator in the development of severe disease and higher levels have been identified in patients with dengue haemorrhagic fever and dengue shock syndrome compared with dengue fever (Hober et al., 1996;Green et al., 1999;Gagnon et al., 2002). TNFα secreted from DENV-infected immune cells leads to fluid leakage by activating vascular endothelial cells (Anderson et al., 1997), it activates the coagulation cascade (Hober et al., 1993) and induces production of vasoactive substances (Anderson et al., 1997).
The effect of DNJ-derived iminosugars on DENV-infected primary DCs has not yet been investigated. In addition, the effect of DNJ-derived iminosugars on TNFα secretion by DENV-infected DCs is yet to be explored. This study examined the antiviral effects of iminosugar treatment on DENV-infected immature dendritic cells (imDCs), as well as the impact on TNFα secretion.

Viruses and cells
DENV2 strain 16681 (a gift from G. Screaton, Oxford, UK) was propagated in the C6/36 insect cell line (a gift from the Armed Forces Research Institute of Medical Sciences, Thailand). Human monocytes were isolated from buffy coats (NHS Blood and Transplant, surplus to clinical requirements) as described previously (Miller et al., 2008) and cultured in RPMI-1640 media (with 10% heat-inactivated foetal bovine serum) for 4 days, supplemented with 50 ng/ml IL-4 and 100 ng/ml granulocyte colony stimulating factor (GM-CSF, Peprotech) to differentiate the cells into imDCs. The use of human blood was approved by the NHS National Research Ethics Service (09/H0606/3). Flow cytometry was performed to confirm the phenotype of imDCs. Cells were scraped from the flasks, washed in cold PBS, re-suspended in FACS block (0.5%w/v BSA, 5 mM EDTA in PBS, 5%v/v heat inactivated human serum) and placed on ice for 30 min and stained with primary antibody for 45 min on ice. 10 μg/ml anti-DC-SIGN (R&D, #MAB161) and 10 μg/ml anti-CD86 (BioLegend, #305401) antibodies were compared to 10 μg/ml IgG2b (R&D, #MAB004) isotype antibody. 50 μg/ml anti-CD14 (R&D, #MAB3832), 20 μg/ml anti-CD16 (BioLegend, #305301), 10 μg/well anti-CD83 (BioLegend, #305307) and 20 μg/ml anti-MR antibodies (clone 3.29.B1; a gift from Luisa Martinez-Pomares, University of Nottingham) were compared to 50 μg/ml IgG1 (clone D1.3) isotype control. Cells were re-suspended in 2 μg/ml secondary antibody (PE-Goat anti-mouse IgG (F(ab')2 (eBioscience)) after two washes with FACS buffer and kept on ice for 45 min. Cells were subjected to two further washes before re-suspending in FACS buffer. Data were collected using a FACS Calibur machine (BD Biosciences) with approximately 10,000 gated cells per population and analysed using BD CellQuest™ software. Histograms were used to determine the expression of the surface receptors by comparing staining with target antibodies with the relevant isotype control.

Cytotoxicity assay
Cytotoxicity of iminosugars was assessed by assaying for metabolic activity using a CellTiter 96® AQueous One Solution Cell Proliferation Assay (Promega) as per the manufacturer's instructions. Samples were incubated for 1-4 h (37 • C, 5% CO 2 ) and the absorbance was measured at 490 nm (A 490 ) on a SpectraMax M5 microplate reader (Molecular Devices). Absorbance was normalized to untreated controls after subtracting the blank reading. CC 50 was defined as concentration of iminosugar at which cell viability was 50%.

Antiviral assays on imDCs
ImDCs were infected with DENV2 strain 16681 at a multiplicity of infection (MOI) of 1 for 90 min at room temperature, then virus containing supernatant was replaced with RPMI-1640 medium (Cambrex) supplemented with 10% foetal bovine serum and serial dilutions of iminosugar. Cells were incubated for 48 h at 37 • C with 5% CO 2 . Then the cell culture supernatant was removed, centrifuged for 5 min at 400×g to pellet any cells, and stored at − 80 • C. Cells were fixed with 4% paraformaldehyde and permeabilised.
Infectious viral titers were obtained by viral plaque assay (Miller et al., 2008) on LLC-MK 2 monkey kidney cells (limit of plaque detection of 33 plaque forming units [PFU]/ml). The percentage of infected cells was measured by immunofluorescence (Miller et al., 2008) Tumour necrosis factor alpha 50% and 90% inhibitory concentrations (IC 50 and IC 90, respectively) were calculated based on experiments with imDCs derived from at least three donors. The selectivity index (SI) was calculated as follows: SI = CC 50 /IC 50 . Total secreted viral genomes were quantified by extracting viral RNA from cell culture supernatant using the Directzol RNA mini prep kit (Zymo Research). A one-step quantitative reverse transcription polymerase chain reaction (qRT-PCR) for DENV NS5 was performed using the Thermo Scientific Verso one-step RT qPCR kit. A 2 μL volume of total template RNA in a final 20 μL reaction was run at 48 • C for 30 min (reverse transcription), 95 • C for 15 min (initial denaturation) followed by 40 cycles of 95 • C for 15 s, 55 • C for 60 s and 60 • C for 25 s before reading the plate. All samples were run in technical duplicates and compared to a standard curve of high titre viral RNA isolated from C6/ 36-grown DENV-2 for quantification. Forward (3 ′ -ACAAGTCGAA-CAACCTGGTCC) and reverse (3 ′ -GCCGCACCATTGGTCTTCTC) primers for dengue NS5 were used. The qRT-PCR was performed using an Applied Biosystems 7500 Fast Real-Time PCR System. Specific infectivity was calculated as follows: Specific infectivity = infectious virus in pfu/total viral genome copies. The value was normalized to the untreated (no iminosugar) control.

Functional TNFα assay
Functional TNFα levels were measured using HEK-Blue™ TNFα reporter cells (Invivogen, catalogue #hkb-tnfdmyd) cultured per manufacturer's instructions. Stimulation of these cells by TNFα results in activation of the NF-κB inducible promoter and production of secreted alkaline phosphatase (SEAP). To detect biologically functional TNFα, reporter cells were stimulated with cell culture supernatants from DENV-infected and/or iminosugar-treated imDCs for 24 h. SEAP secretion was quantified using QUANTI-Blue™ (Invivogen) as per manufacturer's instructions, by measuring absorbance at 645 nm (A 645 ) using a NOVOstar microplate reader (BMG Labtech). Comparison with a standard curve of human recombinant TNFα (Peprotech, catalogue # 300-01A) was used to calculate the functional cytokine levels in the cell culture supernatant.

Detection of free oligosaccharides
Free oligosaccharides (FOS) accumulate due to inhibition of the ER α-glucosidases and quantification of the exact type of FOS provides an indication of αGluI and αGluII inhibition in vitro. Uninfected imDCs or imDCs infected with DENV2 at MOI 1 were treated with 0.3, 3.16 or 31.6 μM concentrations of iminosugars or left untreated. After 48 h postinfection, cells were lysed and analysed by normal phase (NP)-highperformance liquid chromatography (HPLC) for detection of FOS (Miller et al., 2008). For enzymatic determination of individual glucose-terminating and galactose-terminating FOS species, samples were subjected to digestion with β-galactosidase at 50 U/ml at 37 • C overnight to remove the terminal galactose residues. Peak area was used to assess molar quantity in comparison to standards of known identity and quantity. All FOS values were normalized based on the protein concentration.

Statistical analyses
Data were analysed using GraphPad Prism version 9.2.0 (GraphPad Software, Inc). For plaque assay data, the IC 50 was calculated by fitting a logarithmic four-point sigmoidal curve to each data set. One-way ANOVA with Dunnett's multiple comparisons test was used to compare functional TNFα data among iminosugar treatments. Mean values for each set of biological replicates are given in figures with error bars representing standard error of the mean (SEM).

Primary monocyte-derived imDC model
Differentiated imDCs showed increased expression of CD83, CD86, DC-SIGN and mannose receptor (MR) compared to day one monocytes ( Figure A1). Previous literature suggests that imDCs do not support antibody dependent enhancement (ADE) in DENV infection (Marovich et al., 2001;Boonnak et al., 2008). Therefore, we confirmed DENV infection of imDCs in a non-ADE setting. The differentiated imDCs were susceptible to DENV infection as suggested by expression of DENV envelope protein (detected with monoclonal antibody 3H5) 48 h following in vitro infection with the DENV2 strain 16681 (Fig. 1b). Consistent with productive infection, the DENV envelope protein was expressed in 31 ± 5.5% imDCs; secreted virus was detected in the cell culture supernatant by qRT-PCR at 1.4 × 10 8 ± 1.7 × 10 7 genome copies/ml and secreted infectious virus by plaque assay at a titre of 2.4 × 10 5 ± 3.3 × 10 4 pfu/ml.

Iminosugars have minimal cytotoxicity in imDCs
All DNJ-derived iminosugars had low cytotoxicity (Table 1), comparable to that seen in primary macrophages after 48 h of treatment (Fischer et al., 1995). EOO-DNJ treatment did not show significant cytotoxicity even at a concentration of 1 mM. All iminosugars were used within the non-toxic range (at a concentration below the CC 10 ) for further experiments. NN-DGJ was used at concentrations not exceeding 31.6 μM due to variable cytotoxicity at higher concentrations; for example, 100 μM resulted in 82.8 ± 31.67% imDC cell viability after 48 h. All other iminosugars were used at concentrations up to 100 μM.

Iminosugar treatment reduces secretion and specific infectivity of DENV produced in imDCs
To investigate whether iminosugars had an antiviral effect against DENV in imDCs, cells were infected with DENV for 90 min and then treated for 48 h with iminosugars at a serial dose titration from 0.3 to 100 μM. The percentage of infected imDCs reduced with DNJ-based iminosugar treatment in a dose-dependent manner (data not shown).
All DNJ-derived iminosugars tested were antiviral, reducing the amount of infectious virus (Fig. 1c) and total number of DENV virions (Fig. 1d) released. In addition, specific infectivity was significantly reduced at iminosugar concentrations ≥31.6 μM (Fig. 1e). In contrast, treatment with <31.6 μM NN-DGJ did not significantly affect the secretion of DENV from imDCs, nor the specific infectivity of secreted virus.

DNJ-derived iminosugar treatment reduces secretion of functional TNFα levels by DENV-infected imDC
Given the important role of TNFα in DENV pathogenesis, we evaluated the effect of iminosugar treatment on TNFα secretion from DENVinfected imDCs. A functional assay was used to ensure that any TNFα detected was biologically functional, as iminosugar treatment could potentially disrupt the folding of the cytokine, a glycoprotein. DENV-infected imDCs secreted functional TNFα (Fig. 2a) and DNJ-derived iminosugars reduced functional TNFα levels in a dose-dependent manner (Fig. 2b-d). Depending on the iminosugar, TNFα levels were significantly reduced at concentrations of 31.6 μM and/or 100 μM. NN-DGJ treatment did not reduce functional TNFα levels significantly in DENV-infected imDCs (Fig. 2e).

ER αGluI inhibition by DNJ-derived iminosugars increases the antiviral effect and reduces DENV-induced TNFα levels
To better understand the mechanism resulting in the reduction in specific infectivity of DENV and the inhibitory effect of high concentrations of DNJ-derived iminosugars on secreted TNFα levels, iminosugar concentrations at which α-glucosidase inhibition occurs were determined. Accumulation of Glc 1 Man 4 GlcNAc 1 FOS species occurs due to αGluII inhibition and Glc 3 Man 5 GlcNAc 1 species are detected following αGluI inhibition (Alonzi et al., 2008). Measurement of FOS derived from DENV-infected imDCs treated with iminosugars revealed that 2THO-DNJ, EOO-DNJ and NN-DNJ inhibited αGluII at all tested concentrations (Fig. 3a), including 0.3 μM, the concentration corresponding to the onset of reduction in viral secretion. The inhibition of αGluI occurred at 31.6 μM for DNJ-derived iminosugars (Fig. 3b); the concentration at which viral secretion was reduced to less than 20% of untreated samples. Therefore, reductions in both the specific infectivity of secreted virus and functional TNFα levels became apparent at concentrations where αGluI was inhibited. Thus, although αGluII inhibition is sufficient to reduce viral secretion from imDCs, these data suggest that αGluI inhibition could have an additive effect on the antiviral activity and also affect TNFα levels in DENV infection.

Discussion
To our knowledge, this study provides the first evidence of antiviral effects of DNJ-derived iminosugars in DENV-infected primary human monocyte-derived imDCs, a clinically relevant cell type for natural DENV infection (Cheong et al., 2011;Nakano et al., 2009;Schmid and Harris, 2014;Marovich et al., 2001;Ho et al., 2001). In vitro models for antiviral drug evaluation should mimic human infection as closely as possible. Despite limitations of cell culture based experiments, primary cells provide an attractive option to understand cellular responses to infection. Our experiments performed on primary imDCs provide important insights into antiviral and immunomodulatory properties of the host-directed DNJ compounds in DENV infection.
In this study, inhibition of αGluII by DNJ-derived iminosugars was sufficient to induce antiviral effects in imDCs, in keeping with findings for the selective αGluII inhibitor TOP-DNJ in MDMΦs (Kiappes et al., 2018). However, key differences were observed in the antiviral effects elicited by iminosugar concentrations inhibiting only αGluII, compared to concentrations at which both α-glucosidases were inhibited. The reduced specific infectivity demonstrates that virions secreted by imDCs treated with iminosugar concentrations inhibiting αGluI have significantly lower infective potential than those produced in imDCs under αGluII inhibition. Iminosugars were previously thought to have insignificant effects on specific infectivity of secreted DENV (Sayce et al., 2016), contrary to the reduced specific infectivity of human immunodeficiency virus (HIV) produced under NB-DNJ treatment in the absence of a reduced viral output (Warfield et al., 2020). However, our findings in imDCs suggest that DNJ-derived iminosugars have a dual effect, significantly reducing DENV secretion at all antiviral concentrations, and reducing specific infectivity at concentrations higher than 10 μM which inhibit αGluI. Exploitation of αGluI inhibition provides an attractive option for an enhanced antiviral effect and this concept is further supported by the in vivo data available in mice (Miller et al., 2012). Inhibition of αGluI by a single high dose of the iminosugar derivative MON-DNJ (N-(9-methoxynonyl)-1-deoxynojirimycin) prevented death of DENV-infected mice in a lethal antibody-dependent enhanced mouse model even when animals were treated after 48 h post-infection (Miller et al., 2012). These data have important therapeutic implications during dose selection for clinical trials. In addition to the antiviral effects, iminosugar concentrations inhibiting αGluI were associated with reduced TNFα secretion. The TNFα inhibitory effect seen with iminosugar treatment could be due to the reduction of DENV production by imDCs or due to an independent effect on TNFα secretion. A previous study of MON-DNJ revealed that this DNJ-derivative inhibited the inflammatory response elicited by multiple pathogens and suggested that the effect was modulated via the unfolded protein response (Sayce et al., 2021). The mild reduction of TNFα observed with NN-DGJ treatment was not statistically significant or dose-dependent as expected with this compound. The current paper does not attempt to understand the exact mechanism of the reduction in TNFα production observed with iminosugar treatment in imDCs, but experiments designed to answer this question are of interest for the future.

Conclusions
The DNJ-derived iminosugars NN-DNJ, 2THO-DNJ and EOO-DNJ reduce secretion of virus from DENV-infected imDCs in a dosedependent manner. In addition, at concentrations inhibiting αGluI, these compounds reduce the infectivity of the virus that is released, and also reduce TNFα secretion. Taken together, these data suggest that αGluI inhibition by DNJ-derived iminosugars mediates important antiviral effects in DENV-infected imDCs, beyond the reduction in viral secretion achieved with αGluII inhibition alone.

Funding
This work and author JLM was funded by the Oxford Glycobiology endowment. NP was funded by the Commonwealth Scholarship Commission. JB and BET were supported by the Wellcome Trust [grant numbers 203853/Z/16/Z and 105402/Z/14/Z, respectively]. Fig. A1. Phenotyping data for monocyte-derived dendritic cells. Monocytes were isolated from buffy coats and differentiated into day (D) 4 imDCs. D1 monocytes and D4 imDCs were stained using antibodies to detect 6 different surface markers (CD14, CD16, MR, DC-SIGN, CD83 and CD86). Cells were subjected to flow cytometry and the above histograms demonstrate the expression of these surface markers (in blue) compared to the isotype controls (grey filled). The y-axis represents the number of events and the x-axis the fluorescence intensity. These are representative histograms from two independent experiments. MR; mannose receptor, DC-SIGN; dendritic cell-specific ICAM-3 grabbing nonintegrin.

Declaration of competing interest
Authors NP, JB, DA, BET, JLM and NZ declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.