How Complex Should an In Vitro Model Be ? Evaluation of a Complex 3 D Alveolar Model with Transcriptomic Data and Computational Biological Network Models

To more accurately model inhalation toxicity in vitro, we developed a tetra-culture system that combines lung alveolar epithelial cells, endothelial cells, macrophages, and mast cells in a three-dimensional orientation. We characterized the influence of the added complexity using network perturbation analysis and gene expression data. This allowed us to gain insight into the steady-state profile of the assembled, complete three-dimensional model using all four cell types, and of simpler models of one, two, or three cell types. Gene expression data were analyzed using cause-andeffect biological network models, together with a quantitative network-scoring algorithm, to determine the biological impact of co-culturing the various cell types. In the tetra-culture, macrophages appeared to be the largest contributors to overall network perturbations, promoting high basal levels of oxidative stress and inflammation. This finding led to further optimization of the model using rested macrophages, which decreased the basal inflammatory and cell stress status of the co-culture. We compared transcriptional profiles from publicly available datasets of other in vitro models representing the airways and of healthy human lung tissue with those of our model. We found an increasing correlation between airway models and normal human lung tissue as cell types became more physiologically relevant and the complexity of the system increased. This indicates that the combination of multiple lung-relevant cell types in vitro does indeed increase similarity to the physiological counterpart. This is an Open Access article distributed under the terms of the Creative Commons Attribution 4.0 International license (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution and reproduction in any medium, provided the original work is appropriately cited.


Introduction
The past decades have seen a surge in the development of in vitro models driven by 3R efforts, and models available today range in complexity from single cell type monolayers to multi-cell type, three-dimensional (3D) models (Duval et al., 2017;Edmondson et al., 2014).It is generally acknowledged that more complex systems better reflect the in vivo situation (Braakhuis et al., 2015); however, this general assumption should be verified on a case-by-case basis, and the biological contribution derived from the incorporation of each additional cell type should be evaluated carefully.Overall, in vitro cellular models should be as simple pletely resemble the in vivo organization but that is predictive of in vivo outcomes.This leads to the question: How complex should an in vitro model be?
An in vitro alveolar co-culture model combining alveolar type II epithelial cells (A549 cells) (Lieber et al., 1976), macrophage-like cells (phorbol 12-myristate 13-acetate [PMA]-differentiated THP-1 cells) (Tsuchiya et al., 1980;Schwende et al., 1996), mast cells (HMC-1 cells) (Butterfield et al., 1988), and endothelial cells (EA.hy 926 cells) (Suggs et al., 1986) was developed (Klein et al., 2013(Klein et al., , 2017) ) and its relevance was demonstrated upon challenge with a large concentration range of diesel exhaust particulate matter (Klein et al., 2017;Fizeșan et al., 2018).The aim in developing the tetra-culture was to modify and improve an earlier version of such a system (Alfaro-Moreno et al., 2008) so as to assess the biological effects of aerosol exposure of the relevant lung cells organized in a 3D orientation at the air-liquid interface (ALI), i.e., to mimic the in vivo anatomy and functionality of the human alveolar region.To better represent the in vivo anatomy and to improve the structural integrity, the epithelial cells are seeded on the apical part of a Transwell ® membrane, while the endothelial cells are seeded on its basal side.This co-culture is the "backbone" of the tetra-culture system that is further augmented by addition of macrophages and mast cells to the apical compartment.The culture is then lifted to the ALI, allowing the system to be used to examine the alveolar effects of inhaled substances (Klein et al., 2013(Klein et al., , 2017)).Exposure at the ALI is recognized for the significance and relevance of the data that it generates and is expected to be a key strategy in establishing alternative methods in inhalation toxicology and future research in this direction (Lacroix et al., 2018).
The main goal of the present study was to examine the specific contributions of each of the cell types forming the alveolar tetra-culture at the ALI on the transcriptional level and to gain insight into the steady-state profile of the naïve (i.e., untreated) tetra-culture.Cells were cultured in mono-, bi-, tri-, and tetra-culture, and gene expression profiles were generated and analyzed using a computational network approach to determine the impact of the individual cell types on the molecular level.Based on the findings of this investigation, the tetra-culture system was then further refined and reanalyzed.Finally, to understand whether the increased complexity of the tetra-culture system really offers increased physiological relevance, we compared our profiles with publicly available gene expression profiles from other in vitro airway models and from healthy human lung tissues.

Cell cultures and preparation of the tetra-cultures
A549 lung adenocarcinoma cells, EA.hy 926 endothelial cells, and THP-1 monocyte-like cells were obtained from the American Type Culture Collection (Manassas, VA, USA).HMC-1 mast-like cells were provided by J. H. Butterfield, Mayo Clinic cell-cell communication, drug metabolism, gene expression, and protein synthesis, generally increasing their in vivo relevance (Antoni et al., 2015).
A co-culture model is a cell culture system in which two or more different populations of cells are grown together, allowing for some degree of interaction.Co-cultures have become highly relevant for drug research, toxicology, and cellular biology studies, as well as for risk assessment.In fact, they provide a more representative human in vivo-like tissue model and, unlike in vivo animal models, allow for high-throughput testing and in-depth monitoring of effects on cell-cell interactions (Gordon et al., 2015).Co-culture systems could prove to be useful for a more accurate prediction of human hazards than simple cell culture models currently permit, thereby also meeting global needs for non-animal models and more refined approaches in line with the 3R principles.Moreover, the combination of such complex and more physiologically relevant in vitro models with in silico data will be pivotal to safety assessment if animal testing is to be avoided altogether.Before these systems can be used in industrial or regulatory contexts (e.g., in pharmacological research or toxicological studies), however, they must undergo a rigorous process of optimization and verification to ensure functionality and reproducibility.Generally, such validation of alternative methods is performed by following clear guidelines and includes protocol development and internationally performed round-robin exercises organized by the European Union Reference Laboratory for Alternatives to Animal Testing (EURL ECVAM)1 .Ultimately, successful validation and acceptance of EURL ECVAM recommendations may result in acceptance as a test guideline by the Organisation for Economic Co-operation and Development (OECD).Although a number of in vitro models of the (predominantly large) airways are commercially available (Clippinger et al., 2018), to our knowledge, these are not yet validated, nor are there validated alveolar models.With respect to the model described in this work, we have performed intra-and inter-laboratory validation studies, which are highly recommended prior to submission to EURL ECVAM for formal validation.In addition, to further increase its relevance as an alternative in vitro model, we are currently trying to replace all animal-derived products (e.g., cell culture serum, antibodies for immunostaining, trypsin, etc.).
To create human in vivo-like tissue models that are more physiologically relevant, it is important to ensure that all relevant cellular components are represented in the cell culture system and that the basic physiological functions of the organ are represented (Chary et al., 2018).Components such as epithelium, endothelium, and immune cells are all equal contributors to tissue homeostasis and tissue response.Therefore, it is highly important that spatial distribution, relative abundance, and phenotypic characteristics are maintained to the extent that the effects observed are relevant for the in vivo situation.In other words, a model that perfectly resembles the structure of the in vivo tissue but fails to mimic the in vivo physiological functions has limited use, contrary to an in vitro model that does not com-

mRNA expression analysis using GeneChips ®
Target preparation workflow was performed using the Biomek FXP Target Prep Express liquid handling system (Beckman Coulter, Brea, CA, USA).In brief, 100 ng total RNA was reverse-transcribed to cDNA using the Gene Chip ® HT 3′ IVT PLUS kit (Affymetrix, Santa Clara, CA, USA).The cDNA was then labeled, amplified to complementary RNA (cRNA), and fragmented.Successful fragmentation was confirmed using the Agilent 2100 Bioanalyzer, and hybridization cocktails were prepared for each sample according to the manufacturer's instructions.The final cocktails were hybridized to GeneChip ® Human Genome U133 Plus 2.0 Arrays (Affymetrix) at 45°C for 16 h while rotating at 60 rpm.Arrays were then washed and stained on an FS450 DX GeneChip ® Fluidics Station (Affymetrix) and scanned using a 3000 7G GeneChip ® Scanner (Affymetrix) to generate CEL files containing the raw gene expression data.
The raw data CEL files were preprocessed through a pipeline based on multiple Bioconductor packages developed for the R statistical software environment (Huber et al., 2015).Data quality was first controlled by examining log-intensity plots, normalized unscaled standard error plots, relative log expression (RLE) plots, polyA control boxplots, RNA degradation plots, spike-in control boxplots, and pseudo and raw images, all generated using the affyPLM and affy packages (Gautier et al., 2004;Brettschneider et al., 2008).Arrays that fell below a set of thresholds on these quality control checks were excluded from further analysis.Admitted CEL files were then background-corrected, normalized, and summarized using frozen robust microarray analysis (McCall et al., 2010).Background correction and quantile normalization were applied to generate microarray expression values from all arrays using the custom CDF environment HGU133Plus2_Hs_ENTREZG v16.0 (Dai et al., 2005).A moderated statistical method implemented in the limma R package was used to identify significantly regulated genes (Smyth, 2005).Adequate linear models were defined to estimate effects of interest, and raw p-values were generated for each probe set.The Benjamini-Hochberg false discovery rate method was then applied to correct for multiple testing effects (Benjamini and Hochberg, 1995).Adjusted p-values below 0.05 were considered statistically significant.

NPA and BIF
Gene expression profiles were further analyzed in the context of a collection of hierarchically structured network models describing the molecular mechanisms underlying essential biological processes in healthy lung tissues (Boue et al., 2015).The network models can be visualized and downloaded from http:// causalbionet.com.The network perturbation amplitude (NPA) scoring method leverages high-throughput measurements and a priori literature-derived knowledge in the form of network models to characterize activity changes for a broad collection of biological processes at high resolution (Martin et al., 2014).
The NPA itself aims to reduce high-dimensional transcriptomics data by combining the gene expression fold changes into fewer differential backbone values (between a few dozen and 200).By definition of the two-layer structure, no measurements (Rochester, MN, USA).Cell culture and standard tetra-culture assembly were performed as previously described (Klein et al., 2013;Chary et al., 2019).Briefly, endothelial and alveolar epithelial cells were seeded on BD Falcon cell culture inserts (4.2 cm 2 surface area, 1 μm pore size, high-pore density PET membranes for 6-well plates; BD Biosciences, Basel, Switzerland).First, endothelial cells were seeded on inverted Transwell™ inserts.After 4 h, the plate of Transwell™ inserts was turned over, and the epithelial cells were seeded inside the Transwell™ inserts.After three days of incubation, immune cells were seeded on the epithelial cells in the apical part of the system.After a further 4 h, medium from the apical part of the insert was removed to switch to the ALI.Epithelial-endothelial co-cultures are considered the "backbone" of the model.The individual immune cell types were added accordingly to the "backbone" to obtain the different tri-culture variations.
Alternative assembly, including the addition of rested macrophages as characterized by Daigneault and colleagues (2010) and delayed switch to ALI, was performed as follows: To study the effect of adding rested instead of freshly differentiated THP-1 cells, differentiated THP-1 were seeded (2.4 × 10 5 cells/cm 2 ) on top of the A549 cells on Day 4, directly after overnight incubation with PMA (Sigma, purity ≥ 98.0%, stock solution diluted in ultrapure ethanol, working solution at 20 ng/ml in medium; non-rested cells, M0) or after five days of resting in fresh complete medium after removal of PMA (five-day rested cells, M5).HMC-1 cells were added to the apical side of the inserts (1.2 × 10 5 cells/cm 2 ) together with the differentiated THP-1 cells to complete the tetra-culture system.Following attachment of differentiated THP-1 (M0 or M5) and HMC-1 cells (4 h after seeding), the medium was removed from the upper compartment, and the tetra-culture was maintained overnight at the ALI, with 2 ml medium containing 1% fetal bovine serum in the basolateral compartment (ALI0).In parallel, tetra-cultures assembled in the same manner were kept in submerged conditions for an additional 24 h after the seeding of immune cells before switching to the ALI (ALI24).Samples from both groups (ALI0 and ALI24) were collected for RNA extraction after 24 h at the ALI.Three independent biological replicates were generated for each condition.

Sample collection and RNA isolation
Two modes of sample collection were applied: (1) All cell types present were pooled (P), or (2) the endothelial cells were collected from the bottom of the Transwell ® inserts (B, basal collection) and processed separately from the content of the insert (A, apical collection).Cells were collected with a cotton bud into 700 µl QIAzol (QIAGEN, Hilden, Germany).A full randomization was applied across all samples prior to sample processing to minimize batch effects.Total RNA, including microRNA, was isolated using a QIAGEN miRNeasy Mini Kit.The quantity of purified RNA was determined using a NanoDrop™ ND8000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA), while the quality of the RNA was evaluated using an Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA).Only RNA samples with an RNA integrity number ≥ 6 were processed further.the network in extracting the NPA scores (K), significant O and K p-values (typically below 0.05) indicate that the network is specifically perturbed.
The key contributors to the perturbation, termed leading nodes, are by definition the nodes that make up 80% of the NPA score.This score accounts for both the differential backbone values themselves and the centrality of the nodes in the functional layer.
To assess the amplitude of perturbation induced by the different culture conditions defined in this study, the NPA was computed for all contrasts (Tab. 1, 2) using a pre-selection of subnetwork models.NPA scores were then aggregated into a unique value, the biological impact factor (BIF), summarizing the overall perturbation of the biological system.
The mRNA array dataset is available in the Arrays Express repository 2 (ID: E-MTAB-7072), while all leading nodes are reported in the supplementary files 3 .corresponding to nodes in the functional layer are available.The differential backbone values will therefore be determined by a fitting procedure that infers values that best satisfy the directional and signed relationships contained in the backbone model while being constrained by the experimental data (the log 2 fold changes for genes).NPA summarizes into a single number how well and how strongly the inferred values match the network, and the resulting score serves as a basis for further interpretation.As a significant NPA score (with respect to the biological variation) alone does not indicate whether the gene fold changes specifically reflect the network structure itself, so-called companion statistics, O and K, are employed.By assessing the adequacy of the downstream gene (transcript layer) assignment to the nodes of the backbone model by reshuffling the gene labels in the transcript layer (O) and by assessing the importance of the cause-and-effect relationships encoded in the backbone of

Alternative assembly protocol
Three alternative assembly protocols were examined with respect to the steady-state of the assembled tetra-culture.
(1) To study the effect of adding rested instead of freshly activated macrophages (contrast #6), THP-1 cells were washed to remove PMA, supplemented with fresh medium, and returned to the incubator for five days prior to being added to the epithelial-endothelial co-culture with the mast cells.The samples collected for gene expression analysis included those from the original assembly protocol (not rested, termed "M0") and those prepared with the rested macrophages (termed "M5").
(2) To examine the impact of a delayed lift to the ALI (contrast #7), tetra-cultures were prepared following the original protocol (i.e., including a switch to ALI once the immune cell types were added) (termed "ALI0"), or cultures were assembled and returned to the incubator for 24 h before the switch to ALI was performed (termed "ALI24").
(3) Both protocol modifications were combined to produce tetra-cultures with rested macrophages that were switched 24 h after assembly (contrast #8).Three independent biological replicates were generated for each condition.

In vitro and in vivo correlation
The gene expression profiles of human lung tissue in three datasets were used as reference gene expression profiles: − GSE19804 (Genome-wide screening of transcriptional modulation in non-smoking female lung cancer patient in Taiwan, published on NCBI on 30 January 2011)

Contribution of individual cell types
We first investigated the effects of adding the endothelial cells to the epithelial cells and, subsequently, we studied the impact of the immune cells (macrophages and mast cells), adding them either individually or in combination to the epithelial-endothelial co-culture (Fig. 1).The samples collected for RNA extraction/ microarray analysis and the examined contrasts are listed in Table 1.Performing the contrasts following this approach provides an indication of the impact that the addition of each cell type has on epithelial cells (contrast #1), endothelial cells (contrast #2), and the whole culture (contrasts #3, #4, and #5).It is worth noting here that apical and basolateral cells were harvested separately to specifically address the effect of adding either cell type.
Considering that a single cell type is harvested and analyzed, this approach enabled the evaluation of mutual effects in the context of the epithelial-endothelial co-culture.When immune cells were added, however, the effects at the single cell-type level could no longer be analyzed, so cells harvested from both the apical and basal sides were pooled to enable the evaluation of the effects on the whole culture.NPA scores and relative BIFs (RBIF), when possible, were determined using a selection of biological network models that capture the biology of cell stress (CST) responses, cell proliferation (CPR), cell fate (CFA) (i.e., apoptosis, necroptosis, autophagy, and senescence), inflammation (IPN), and tissue repair and angiogenesis (TRA).protein families were inferred to be upregulated, whereas the molecules known to halt cell cycle progression (MDM2, TP53, RB1, and cyclin-dependent kinase inhibitors) were inferred to be downregulated in the epithelial compartment.The directionality of the inferred impact on the cell cycle regulators was opposite in the endothelial compartment, indicating promotion of senescence rather that cell cycle progression (Bandara et al., 1994;Duronio and Xiong, 2013).The leading node analysis of the cell cycle network model, which was not impacted in the endothelial compartment, further supported cell cycle progression in the epithelial compartment (Fig. 4).G1/S transition of the cell cycle was inferred to be increased, resulting from activation of E2F proteins, and cell senescence was inferred to be decreased.Probes in the raw CEL files were summarized using the revised Entrez-based probe annotation (Dai et al., 2005) and normalized using frozen robust microarray analysis (McCall et al., 2010).Quality controls, including log-intensities, normalized unscaled standard error, RLE, median absolute value of RLE, and pseudo-images, as well as raw image plots, were computed with the affyPLM package (Bioconductor suite) (Bolstad et al., 2005).Gene expression profiles were averaged per experimental group and dataset.Pearson correlations of gene expression profiles between in vitro and in vivo gene expression profiles were estimated.A bar plot of the correlation coefficients was generated using the R package ggplot2 (Wickham, 2016).

Epithelial-endothelial co-culture
First, we sought to evaluate the effects of co-culturing epithelial and endothelial cells on both cell types.Epithelial cells and endothelial cells were harvested independently from the co-culture and then compared with the corresponding mono-culture using a set of computational networks representative of CST, CPR, and CFA biology.
Only a few of the selected networks were impacted when epithelial and endothelial cells from co-cultures were compared with their respective mono-cultures (Fig. 2).The leading node analysis of the senescence network model indicated that the key contributors to the observed network perturbations were impacted differently in the epithelial and endothelial compartments in co-culture conditions (Fig. 3).In particular, cell cycle-promoting molecules, such as TFDP1 and members of the E2F and CDK The leading node analysis of the mast cell activation network model suggested upregulation of pathways that lead to mast cell activation (Fig. 6).Inferred TLR2 activation leading to the activation of the PLA2 family of phospholipases was specific to cultures containing macrophages, whereas PLCG1 (phospholipase C) and the PRKC family were among the leading nodes only in cultures containing mast cells.Increased production of various leukotrienes were common leading nodes for both immune cell types, lending support to the notion of completed monocyte-to-macrophage differentiation and acquisition of macrophage characteristics (Coffey et al., 1994;Covin et al., 1998).
Not surprisingly, we noted a distinct impact of adding macrophages to the epithelial-endothelial co-culture on processes captured in the macrophage signaling network.The leading node analysis suggested upregulation of CD36 and ALOX5, leading to increased leukotriene production and phagocytosis, respectively.TNFα, TLR2, and NF-κB signaling were also inferred to be increased (Supplementary File 2 3 ).Together, these mechanisms suggest a polarization of macrophages toward the

Impact of the presence of macrophage-like and mast cells on epithelial-endothelial co-culture
To evaluate the impact of the presence of immune cells, the previous investigation was extended to tri-cultures, composed of the epithelial/endothelial co-culture with the addition of either macrophages or mast cells, and the tetra-culture.The analysis also included additional relevant biological networks, such as the IPN, TRA, and vascular inflammatory processes (VIP) networks.
The highest biological impact was observed when all four cell types were cultured together (RBIF = 100%, reference; Fig. 5A).With very few exceptions, almost all subnetworks related to CPR, cell senescence and inflammatory processes were strongly perturbed when the complete tetra-culture was considered (Fig. 5B).This pattern of network perturbations was closely mirrored by results from both tri-cultures.The addition of macrophages, however, appeared to perturb the networks of interest more strongly than the addition of mast cells, i.e., the network perturbations elicited by adding only the mast cells were fewer in number, with overall lower NPA scores.

Evaluation of alternative assembly protocol
Taken together, the data support the notion that the tetra-culture is a stable, homeostatic system at the time of or shortly following assembly.There was a strong indication, however, of active cell-cell interaction through the release of soluble mediators that continued to (mostly) positively affect cell growth in the complete tetra-culture.Moreover, the addition of macrophages to epithelial/endothelial co-cultures appeared to elicit a pro-inflammatory phenotype.At the same time, the data indicated a high basal level of oxidative stress, exemplified by perturbations in the NFE2L2 signaling, oxidative stress, and macrophage activation networks.These findings led to the hypothesis that the tetra-culture would respond to respiratory irritants and inhalation toxicants that are known to activate inflammatory pathways and/or induce ROS and ROS-mediated signaling.However, in the case of small to moderate effects, it could be difficult to discern the biological effects from the baseline levels.
Therefore, we decided to assess three modifications to the original assembly protocol (Klein et al., 2013).First, introducing a five-day resting period prior to adding PMA-activated THP-1 cells (termed M5) to the co-culture was thought to enhance THP-1 differentiation to macrophages while lowering basal oxidative stress levels.Second, we hypothesized that a 24-hour delay between tetra-culture assembly and the switch to ALI (termed (classically activated) M1 phenotype, capable of promoting a pro-inflammatory environment by producing various cytokines and inducing tissue damage through the generation of nitric oxide and reactive oxygen species (ROS) (Wang et al., 2014;Verreck et al., 2004).Inferred increases in transcriptional activity of the NF-κB complex, TNFα receptor, and MAPK in the epithelial innate immune response network point toward an acute inflammatory response, which could potentially threaten the stability of the co-culture by increasing epithelial and/or endothelial permeability (Dreymueller et al., 2012;Arndt et al., 2011).Of note, key drivers of endothelial cell activation subnetwork perturbations partially overlapped with those indicative of epithelial cell-mediated pro-inflammatory signaling and further support the hypothesis of an epithelial cell/macrophage-mediated inflammatory response (Supplementary File 2 3 ).In comparison, the addition of mast cells, although inducing a qualitatively similar network response, did not exert such a strong impact on the epithelial/ endothelial co-culture.
Minor perturbations were noted in the CST subnetworks following the addition of immune cells to the epithelial/endothelial co-culture.Overall, the impact of macrophages appeared to be greater than the impact elicited by mast cells, and once again, the predicted mechanisms underlying CST responses were similar for both immune cell types (Supplementary File 2 3 ).ALI24) may allow for a better acclimation of the four cell types to the co-culturing conditions and potentially reduce stress signaling.Finally, both modifications were combined to generate the tetra-culture, and all modifications were examined with respect to their impact on the culture system (Tab.2).
The highest biological impact was obtained when the cultures were lifted to the ALI 24 h after the addition of the unrested macrophages (ALI24) (Fig. 7A).This was partially mirrored when rested macrophages were added simultaneously (M5-ALI24).In contrast, adding rested macrophages (M5) not only appeared to exert a lower biological impact, but was also characterized by a considerably different underlying biology compared with the reference contrast (ALI24).
The NPA heatmap confirmed that the delayed switch to the ALI had a stronger impact on the culture system than the addition of rested macrophages, and this effect appeared to predominate when both modifications were applied simultaneously (Fig. 7B).
The mechanisms related to inflammatory processes were the main contributors to the overall biological impact.Scoring of the macrophage signaling network model with transcriptomic data and the leading node analysis indicated that MIR146a was upregulated when macrophages were allowed to rest prior to culture assembly.In the network model, this feeds into the inhibition of TRAF6 and IRAK1 signaling.Similarly, TLR2 and IRAK4, which also activate IRAK1 and TRAF6 and eventually the NF-κB pathway, were inferred to be downregulated in the M5 cultures compared with the original culture setup.In contrast, the essential molecules leading to and within the NF-κB pathway were inferred to be upregulated following delayed ALI and, to a lesser extent, when delayed ALI cultures contained rested macrophages (Fig. 8).Consequently, IL-1β and TNFα signaling were inferred to be downregulated in the M5 cultures and upregulated in the ALI24 and ALI24-M5 cultures compared with the original culture setup.The network analysis also suggested downregulation of ALOX 5 and leukotriene B4 abundance after macrophage resting (Fig. 8).Further details of the affected signaling pathways within the macrophage signaling network model are illustrated in Figure 8.
To summarize, the directionality of the majority of the common leading nodes in the macrophage signaling network model were inferred to be opposite in the cultures with rested macrophages and delayed switch to the ALI compared with the original assembly protocol.In essence, the addition of rested macrophages resulted in reduced activation, whereas delayed ALI stimulated these signaling pathways.The ALI24 effect was partially mirrored in ALI24 cultures with rested macrophages.
Macrophage resting also elicited minor but significant perturbations in the CST networks.Perturbations in the oxidative stress network were driven mostly by an inferred increase in catalase activity and decreased MAPK1 and PRKCA signaling.MAPK Erk1/2, NF-κB, and FOS/JUN signaling were also inferred to be decreased.Based on these findings, the overall oxidative status of the culture appeared to be decreased (Supplementary File 3 3 ).
Maintaining the complete tetra-culture in submerged conditions for 24 h after the immune cells were added had a strong impact on the cellular stress level of the culture, as oxidative, Fig. 5: Biological network analysis of the effects of the immune cell addition A: RBIF for immune cell contrasts (see Tab. 1).For each comparison, the δ value (−1 to 1) indicates how similar the underlying network perturbations are with respect to the reference (REF, RBIF = 100%).A value of 1 indicates that all networks are perturbed by the same mechanisms.B: Heatmap of NPA scores for immune cell contrasts summarizing subnetworks that were significantly perturbed in at least one contrast.A network is considered perturbed if, in addition to the significance of the NPA score with respect to the experimental variation, the two companion statistics (O and K), derived to inform on the specificity of the NPA score with respect to the biology described in the network, are significant.Significantly perturbed networks are indicated by * (i.e., O and K statistic p-values < 0.05 and significant with respect to the experimental variation).CFA, cell fate; CPR, cell proliferation; CST, cell stress; IPN, inflammatory process network; TRA, tissue repair and angiogenesis; VIP, vascular inflammatory processes system.We also considered available data from primary cell-derived 3D cultures and other 2D cultures.
When comparing gene expression profiles of 3D organotypic models of different regions of the respiratory tract (e.g., nasal) with profiles in the human lung, the degree of similarity was low (Fig. 9).Surprisingly, on the transcriptional level, simple mono-cultures exhibited a higher similarity to the in vivo tissue than certain 3D models.Importantly, however, the analysis demonstrated that combining different cell types, such as epithelial, endothelial, and immune cells, results in transcriptional profiles that strongly correlate with those reported for normal human lung tissues (Fig. 9), with tri-cultures and tetra-cultures showing the highest level of similarity to normal human lung tissue.The reason why this correlation surpasses even that between 3D organotypic models and healthy human lung tissues could be explained by the lack of endothelial and/or immune cells in the 3D buccal, gingival, nasal, and bronchial tissues assessed.It is worth noting that tri-and tetra-culture gene expression profiles were both similar to profiles in human lung tissue, raising the question of whether the addition of a not fully functional mast cell line (HCM-1 cells do not express FcεRI (Nilsson et al., 1994)) is truly necessary, as the degree of similarity does not increase remarkably, whereas the complexity of the model does.
osmotic, and hypoxic stress responses were activated.Increased NADPH oxidase complex, which could potentially increase the levels of ROS, was predicted.Similarly, the catalytic activity of G6PD, another antioxidant enzyme, was predicted to decrease.Lower levels of this enzyme, which participates in the pentose phosphate pathway, may have decreased the levels of the co-enzyme NADPH, leading to glutathione depletion, as described previously (Thomas et al., 1991).
Together, these results suggest that the 24-hour rest in submerged conditions prior to the switch to ALI exerted a detrimental effect on the stability of the system and enhanced the pro-inflammatory status of the tetra-culture.This could be partially attenuated by the use of rested macrophages.

In vitro to in vivo comparison
To further evaluate and benchmark the relevance of the observed cell-cell interactions and behavior of the tetra-culture model, we performed a comparison at the transcriptional level with publicly available datasets of conventional lung in vitro models and of healthy human lung tissues.We included the mono-cultures that are used as the starting point to assemble the complete tetra-culture and the bi-, tri-, and tetra-cultures examined in this evaluation to study the effects of adding each cell type on the complete

Discussion
Rodents are widely used as in vivo models in respiratory toxicology.There are, however, key differences in the anatomy of the pulmonary system of humans and mice or rats.These simple anatomical differences in size and diameter of the respiratory tract can lead to meaningful variance in particle deposition, questioning the relevance of rodents as suitable models to examine the toxicity of particulates (Warheit and Hartsky, 1990;Schlesinger, 1985).Furthermore, differences between the human and rodent immune systems may question the relevance of studying pathophysiological reactions to inhalation exposures in rodent models (Mestas and Hughes, 2004).The time and cost of animal experimentation and the desire to replace animal research with alternative methods explains the drive behind the development of human cell-based models that exhibit the characteristics of their physiological counterparts in both architecture and function.Such models may prove to be valuable tools in the study of respiratory toxicants.Thus, achieving close resemblance of an in vitro system to the tissue it mimics while also capturing the tissue's ability to respond to an exposure stressor requires a certain complexity.However, if the model is too complex, it may become impractical as an alternative test system in an industrial setting.In brief, in vitro models should be as complex as needed and as simple as possible with respect to the number of cell types, the assembly mode, and the number of readouts.
The original model that formed the starting point for the current investigation is complex, but its suitability for studying respiratory toxicants that activate inflammatory pathways has been shown (Klein et al., 2013;Brinchmann et al., 2018).In this study, we attempted to delineate the molecular features of cell-cell crosstalk upon co-culture of up to four cell types that we consider representative of the human alveolar region: alveolar type II epithelial cells, endothelial cells, macrophages, and mast cells.Gene expression data were analyzed within biological cause-and-effect networks using computational tools to gain insight into the steady-state profile of the assembled complete 3D model using all four cell types in comparison to the less complex versions (mono-, bi-, and tri-cultures).Although the culture assembly takes multiple days and requires many manipulations, our data indicate that the tetra-culture is a stable, homeostatic system once assembled.Unexpectedly, the addition of macrophages appeared to influence the co-culture more than the addition of mast cells.
Indeed, in the two previous studies (Klein et al., 2013;Brinchmann et al., 2018), where the tetra-culture was exposed to various concentrations of diesel exhaust particulate matter, translocation of the transcription factor Nrf2 and significantly altered expression of CYP1A1 mRNA in the non-directly exposed endothelial cells were observed (Klein et al., 2017).These observations suggest some exchange between the different compartments of the system, as cell-to-cell communication may happen at the level of the pores.These results are in agreement with the activation of multiple metabolic pathways after the addition of the macrophages.However, the detailed mechanisms of the "cellular crosstalk" that controls the activation of different pathways in the tetra-culture are still under investigation.However, the baseline status of the model was never compared with other models and tissues.Particularly the need for mast cells, which add to the system's complexity, was never demonstrated against simpler versions comprising fewer cell types.HMC-1 cells are deficient in IgE receptor-mediated signaling, a hallmark of allergic hypersensitivity reactions (Sibilano et al., 2014), therefore the relevance of this cell type within the system is questionable.Also, allergic hypersensitivity is only relevant for specific testing questions and not, for example, for assessing the toxic effects of tobacco product aerosols, which predominantly activate macrophages in vivo (Stämpfli and Anderson, 2009).Therefore, the additional complexity brought about by adding mast cells to the system may only be warranted for respective research questions.
Our results (Fig. 10) suggested that the addition of macrophages to epithelial/endothelial co-cultures elicited a pro-inflammatory and oxidative stress phenotype.While such stress responses can be expected when macrophages are co-cultured with epithelial or endothelial cells (Che et al., 2017;Dehai et al., 2014;Loret et al., 2016;Kasper et al., 2017), other studies did not record increases in the expression of genes associated with inflammatory responses or secreted inflammatory mediators (Bengalli et al., 2017;Müller et al., 2009).We reasoned that this baseline activation may ultimately lead to masking of small exposure effects in a test setting, particularly when potential respiratory irritants or toxicants act by activating inflammatory pathways and/or ROS-mediated signaling.Therefore, we investigated modifications to the culture assembly protocol to improve its stability and sensitivity.The use of macrophages that had rested for five days following differentiation suggested a dampening of the previously observed pro-inflammatory phenotype.
In contrast, delaying the switch to the ALI for 24 h seemed to further enhance oxidative, osmotic, and hypoxic stress, suggesting that this modification exerts a detrimental effect on the stability of the tetra-culture.Combining the two modifications (i.e., adding rested macrophages and delaying the switch to the ALI for 24 h) yielded no obvious benefit.In fact, the use of rested macrophages only partially attenuated the heightened stress responses seen when raising the assembled tetra-culture to the ALI only after the 24-h delay.Therefore, future testing approaches for our research questions will involve tri-cultures with rested macrophages.
Finally, we sought to examine how the tetra-culture compares with the in vivo counterpart it mimics, i.e., the human lung.We utilized multiple publicly available gene expression profiles of healthy human lungs for comparison with those of the mono-, bi-, tri-, and tetra-cultures.In addition, we compared gene expression profiles of our culture system with expression profiles of other 3D airway models.Correlation analysis showed that the choice of cell types and the 3D organization are important factors that influence the physiological relevance of a test system.Clearly, the addition of an immune component increases the complexity of the model and may even threaten the stability of the test system, depending on the immune cell type, its activation status, and the sequence of model assembly.Careful evaluation of the consequences of cell-cell interactions in complex cultures,  however, can guide further modifications that yield a stable test system that closely resembles the organ it aims to mimic.
The data presented here speak to the molecular features of the tetra-culture rather than to its functionality or ability to mount adequate biological responses to stressors.Additional tests and functional read-outs should be performed to confirm that the model is fit for its purpose (i.e., the testing of potentially toxic effects of inhaled substances) and will further guide selection of the number and type of cell lines to be included in the model depending on its intended use.

Fig. 1 :
Fig. 1: Tetra-culture composition and data analysis strategyGraphical representation of the different set-ups of the 3D model from the single cultures (either epithelial cells or endothelial cells), to bi-cultures (epithelial cells and endothelial cells), tri-cultures (co-cultures of epithelial and endothelial cells with the addition of either macrophages or mast cells), or the complete tetra-culture (epithelial cells, endothelial cells, macrophages, and mast cells).The different set-ups were used for the NPA analysis.

Fig. 2 :
Fig. 2: Biological network analysis of epithelial-endothelial cell co-culturesHeatmap of NPA scores summarizing subnetworks that were significantly perturbed in at least one contrast (see Tab. 1).A network is considered perturbed if, in addition to the significance of the NPA score with respect to the experimental variation, the two companion statistics (O and K), derived to inform on the specificity of the NPA score with respect to the biology described in the network, are significant.Significantly perturbed networks are indicated by * (i.e., O and K statistic p-values < 0.05 and significant with respect to the experimental variation).CFA, cell fate; CPR, cell proliferation; CST, cell stress

Fig
Fig. 3: The impact of epithelialendothelial co-culture on the senescence network model compared with the respective mono-cultures The graphical representation shows an extract of the senescence network model based on the leading nodes for the comparison of epithelialendothelial co-cultures with the respective mono-cultures.The most contributing leading nodes were prioritized and connected with other leading nodes from the analysis.Some redundant and feedback edges were removed for clarity.The backbone NPA values with directionalities of inferred regulation are shown as bar graphs for each node.Orange/red bars indicate inferred upregulation, and blue bars indicate inferred downregulation.Sample order in the bar graphs: 1) endothelial cells, 2) epithelial cells.The green asterisk indicates that the node is a leading node.Biological expression language (BEL) 4 is used.

Fig. 4 :
Fig. 4: The impact of epithelialendothelial co-culture on the cell cycle network model compared with the respective mono-cultures The graphical representation shows an extract of the cell cycle network model based on the leading nodes for the comparison of epithelial-endothelial co-cultures with the mono-cultures.The most contributing leading nodes were prioritized and connected with other leading nodes from the analysis.Some redundant and feedback edges were removed for clarity.The backbone NPA values with directionalities of inferred regulation are shown as bar graphs for each node.Orange/red bars indicate inferred upregulation, and blue bars indicate inferred downregulation.The green asterisk indicates that the node is a leading node.Biological expression language (BEL) 4 is used.

Fig. 6 :
Fig. 6: The impact of macrophages and mast cells on the mast cell activation network model The graphical representation shows an extract of the mast cell activation network model based on the leading nodes in the different assembly protocols.The most contributing leading nodes were prioritized and connected with other leading nodes from the analysis.Some redundant and feedback edges were removed for clarity.The backbone NPA values with directionalities of inferred regulation are shown as bar graphs for each node.Orange/red bars indicate inferred upregulation, and blue bars indicate inferred downregulation.Sample order in the bar graphs: 1) addition of the macrophages, 2) addition of the mast cells, 3) addition of the macrophages and mast cells.The green asterisk indicates that the node is a leading node.Biological expression language (BEL) 4 is used.

Fig. 7 :
Fig. 7: Biological network analysis of the effect of the modified assembly protocolsA: RBIF for rested macrophages and delayed switch to ALI contrasts (see Tab. 2).For each treatment comparison, the δ value (−1 to 1) indicates how similar the underlying network perturbations are with respect to the reference (REF, RBIF = 100%).A value of 1 indicates that all the networks are perturbed by the same mechanisms.B: Heatmap of NPA for rested macrophages and delayed switch to ALI contrasts summarizing subnetworks that were significantly perturbed in at least one contrast.A network is considered perturbed if, in addition to the significance of the NPA score with respect to the experimental variation, the two companion statistics (O and K), derived to inform on the specificity of the NPA score with respect to the biology described in the network, are significant.Significantly perturbed networks are indicated by * (i.e., O and K statistic p-values < 0.05 and significant with respect to the experimental variation).CFA, cell fate; CPR, cell proliferation; CST, cell stress; IPN, inflammatory process network; TRA, tissue repair and angiogenesis; VIP, vascular inflammatory processes

Fig. 8 :
Fig. 8: The impact of the alternative assembly protocols on the macrophage signaling network model The graphical representation shows an extract of the macrophage signaling network model based on the leading nodes in the different assembly protocols.The most contributing leading nodes were prioritized and connected with other leading nodes from the analysis.Some redundant and feedback edges were removed for clarity.The backbone NPA values with directionalities of inferred regulation are shown as bar graphs for each node.Orange/red bars indicate inferred upregulation, and blue bars indicate inferred downregulation.Sample order in the bar graphs: 1) rested macrophages, 2) delayed ALI, 3) rested macrophages in delayed ALI.The green asterisk indicates that the node is a leading node.Biological expression language (BEL) 4 is used.

Fig. 9 :
Fig. 9: Correlation between normal human lung tissue and selected in vitro 2D and 3D cell culture systems based on transcriptomic profiles Examples of different cellular model systems (x-axis), including both 2D (shades of purple) and 3D (shades of blue) models, are shown.The analysis includes primary (indicated by *) and immortalized cellular models.

2: List of contrasts used to evaluate the effects of five days of rest following induced macrophage maturation and 24-hour delay switching to the ALI
2 https://www.ebi.ac.uk/arrayexpress/ 3 doi:10.14573/altex.1811221sTab.

Tab. 1: Contrasts used to evaluate the effects of different cell type combinations on epithelial cells (apical harvest), endothelial cells (basal harvest), and co-culture (pooled harvest) Contrast Harvest Cellular sample Cellular sample Analysis
Cell types in bold font represent the actual RNA sample.Mac, macrophages; Mas, mast cells