Sex as a biological variable. Our study examined male and female animals, and similar findings are reported for both sexes.

Mouse strains and genotyping. Mice were maintained in ventilated microisolator cages and were fed a standard diet. Mice were provided igloos and nestlets as enrichment items. The LSL-KrasG12D, Prox1-CreERT2, Flt4-CreERT2, and mT/mG strains were described previously (32–35). Established protocols were followed to genotype the mice (17).

Tamoxifen

Newborn mice. We dissolved tamoxifen (25 mg; MilliporeSigma, T5648) in a mixture of ethanol (100 μL; MilliporeSigma, E7023) and peanut oil (900 μL; MilliporeSigma, P2144). To induce Cre-mediated recombination, we fed newborn mice 2 μL of tamoxifen on P0, P1, and P2 using a P20 pipette. Alternatively, tamoxifen was dissolved in sunflower oil and ethanol for a final stock of 20 mg/mL with 5% ethanol. Mice received a subcutaneous injection (100 μg) in the nape of the neck on P0 and P2.

Young adult mice. We dissolved 20 mg tamoxifen in a mixture of 100 μL ethanol and 900 μL sunflower oil (MilliporeSigma, W530285). To induce Cre-mediated recombination, we injected mice (i.p.) with 100 μL of tamoxifen on P22, P24, P26, P28, and P30.

Single-cell isolation, sorting, and sequencing

Lungs from male and female mice were pooled together to isolate LECs from groups of control and KrasG12D mice. The lungs were perfused with heparin/DPBS pH 7.2 to remove red blood cells. The lobes of the lung were isolated and cut into pieces with a scalpel before placing them into C-tubes (Miltenyi Biotec, 130-093-237) containing 5 mL of a lung dissociation cocktail (Miltenyi Biotec, 130-095-927). The C-tubes containing the tissue/enzyme mix were then placed into the sleeve of a gentleMACS Dissociator (Miltenyi Biotec, 130-093-235), which is a benchtop instrument used to dissociate tissues into single-cell suspensions. We ran the preset m_lung_01 program for 8 seconds to dissociate the lung tissue. The C-tubes were then placed in a 37°C incubator under continuous agitation for 45 minutes. Following incubation, the C-tubes containing the tissue/enzyme mix were placed into the gentleMACS Dissociator, and the preset m_lung_02 program was run for 38 seconds to dissociate the cells. The C-tubes were then centrifuged at 300g to collect the sample at the bottom of the tube. The sample was resuspended with 1× Buffer S and applied to a prewet 100-μm cell strainer in a 50 mL conical tube. This process was repeated with a 70 μm and 40 μm cell strainer to ensure that all cells were in a single-cell suspension. The cell suspension was transferred to a 15 mL conical tube and centrifuged at 300g to collect a pellet. The pellet was resuspended in 3% FBS in 1× DPBS for subsequent FACS analysis. The FACSAria II (BD BioSciences) was used for sorting LECs. Cells were sorted at a temperature of 4°C with a 100 μm nozzle at a pressure of 100 psi. Events were first gated on FSC-A versus SSC-A. Then, events were gated on FSC-A versus FSC-W and SSC-W versus SSC-A for gating singles cells. Dead cells were then labeled with DRAQ7 far-red viability dye (Novus Biologicals, NBP2-81126) and events were gated for live cell populations. Finally, cells were gated on GFP-A versus tdTomato-A to sort all GFP-positive LECs. The sorted LECs were then submitted for sequencing. The cells were sequenced by the University of Texas Southwestern Medical Center genomics core. The Chromium Next GEM Single Cell Reagent Kits v3.1 Dual Index (10× Genomics, PN-1000269) was used for library preparation according to the manufacturer’s protocol. Single Cell 3′ v3 Chemistry was used for library preparation. The libraries were sequenced using the NovaSeq 6000 sequencing system (Illumina). Demultiplexing was performed using the bcl2fastq v2.20.0 software (Illumina) and Cell Ranger v7.0.0 (10× Genomics) was used for barcode processing, gene counting, and aggregation.

scRNA-Seq data processing

The Seurat package (v4.4.0) (36) was used to analyze the scRNA-Seq data in R-studio (R v4.3.1). Briefly, raw data from the 10× Genomics Cell Ranger pipeline was loaded into Seurat to generate a UMI count matrix. The count matrix was utilized to create a Seurat object to use in subsequent analysis. While creating the Seurat object, we performed an initial filtration step to remove very-low-quality cells and speed up downstream computational processing. We set the min.cells and min.features arguments to include features that are expressed in at least 3 cells and cells expressing at least 200 genes. After creating the Seurat object, we visualized the QC metrics and performed a second filtration step. Cells that had more than 500 unique features and a percentage mitochondrial count less than 5% were retained. The default parameters were used for control and KrasG12D datasets to normalize the counts, find variable features, scale the data, and perform PCA dimensionality reduction using the NormalizeData(), FindVariableFeatures(), ScaleData(), and RunPCA() functions of Seurat, respectively. To cluster control and KrasG12D datasets, we used Seurat’s FindNeighbors() and FindClusters() functions using 19 PCs and a resolution of 0.3. Dimensionality reduction was performed using the UMAP method implemented in Seurat’s RunUMAP() function with the number of PCs set to 19. The UMAPs were visualized using Seurat’s DimPlot() function. After visualization of the UMAPs, we identified LEC clusters that expressed lymphatic endothelial markers Prox1 and Flt4; blood endothelial cell clusters that expressed the blood endothelial markers Nrp1 and Flt1; and immune cell clusters that expressed high levels of the immune cell marker Ptprc. We also observed additional clusters (in control and KrasG12D datasets) that did not express the lymphatic marker Prox1 uniformly (i.e., low percentage expression) and were composed of cells that on average had fewer numbers of genes (nFeature_RNA) and molecules (nCount_RNA) compared with the other identified clusters. These clusters may represent damaged or senescent cells and were excluded from further analysis. For downstream analysis, we only subset the LEC clusters. Additionally, we subset cells with expression levels of Vegfr3 greater than zero (Flt4 > 0) and levels of Vegfr1 equal to zero (Flt1 = 0). These filtering steps were included to ensure optimal quality and purity of LECs for further downstream analysis. After each subset step, we repeated all the steps from normalization to UMAP dimensionality reduction on control and KrasG12D datasets, only adjusting the parameters for FindNeighbors(), FindClusters(), and RunUMAP(). The clusters were annotated using known LEC subtype–specific markers (21, 22).

Integration of control and KrasG12D datasets

After loading the individual datasets and creating Seurat objects, the control and KrasG12D datasets were merged using the merge() function. We then performed an initial data cleaning step to filter cells based on QC metrics as described in the previous section. The control and KrasG12D datasets were batch corrected and integrated using Seurat’s data integration workflow, which involves performing canonical correlation analysis (CCA) and identifying mutual nearest neighbors (37). After integration, the Seurat workflow involving ScaleData(), RunPCA(), RunUMAP(), FindNeighbors(), and FindClusters() was performed. The clusters were visualized using the DimPlot() function. We then subset the integrated dataset to include only those clusters that expressed Prox1 uniformly (i.e., high percentage expression) and did not contain cells with very low numbers of genes (nFeature_RNA) or low numbers of molecules (nCount_RNA). Additionally, we subset cells with expression levels of Vegfr3 greater than zero (Flt4 > 0) and levels of Vegfr1 equal to zero (Flt1 = 0). After each subset step, we re-integrated, scaled, clustered, and visualized the cells as mentioned above. The clusters were annotated using known LEC subtype–specific markers as described above in scRNA-Seq data processing.

Integration of control and developmental datasets

The control and 4 developmental datasets (27) were loaded, and Seurat objects were created for each dataset. The datasets were merged, filtered based on QC metrics, integrated, and scaled and processed in a similar way as the control and KrasG12D datasets. The clusters were annotated using known LEC subtype–specific markers as described above.

DGE analysis

DGE analysis was performed using the MAST R/Bioconductor package (38). The FindAllMarkers() function in Seurat was used for determining cell type–specific markers for heatmaps, with the min.pct argument set to 0.25, logfc.threshold argument set to 0.25, and test.use argument set to “MAST.” To create a list of all differentially expressed genes in the control and integrated datasets, we used the FindAllMarkers() and FindMarkers() functions, respectively, with the min.pct argument set to 0.1, logfc.threshold argument set to 0.25, only.pos argument set to FALSE, and test.use argument set to “MAST.” The FindMarkers() function was used to obtain the differentially expressed genes for each KrasG12D versus control LEC subtype.

Trajectory inference analysis

The monocle method was used to conduct trajectory inference analysis on control and KrasG12D datasets using the Monocle 3 R package (v1.3.4) (39–41). We used the UMAP coordinates, and cluster labels previously generated by Seurat by converting the Seurat object to the CDS object.

GO analysis

GO analysis was performed using the Metascape web resource (https://metascape.org/gp/index.html#/main/step1) (42). To determine the GO terms from a list of input genes, Metascape employs the hypergeometric test and Benjamini-Hochberg P value correction algorithm (42). To determine the GO terms associated with the control LEC subtypes, we provided Metascape with a list of significantly upregulated genes from each subtype. For analyzing the integrated dataset, we used a list of significantly upregulated or downregulated genes in KrasG12D LECs as an input for Metascape. We plotted the biologically relevant GO terms using GraphPad Prism statistical analysis software (v9.5.1).

Capillary and collecting LEC gene signatures

To determine the capillary and collecting vessel signature scores, we used the AddModuleScore() function in Seurat. The marker genes used to calculate the capillary signature scores were Lyve1, Gjc2, Fndc1, Piezo2, Ccl21a, and Ackr2. The marker genes used to calculate the collecting signature scores were Bgn, Foxp2, Plk2, Sertad1, Fabp4, Ccn1, Cdkn1a, Adamts1, Nos3, Apoe, Ccdc3, Gpx1, Lrg1, and Sgk1. The capillary and collecting signature scores were graphed as violin plots using Seurat’s VlnPlot() function.

Whole-mount immunofluorescent staining

We fixed ears overnight with 1% paraformaldehyde (PFA). Next, we washed the fixed samples with PBS, permeabilized them with PBS plus 0.3% Triton X-100 (PBST), and blocked them overnight with PBST plus 3% donkey serum. After that, we incubated the samples overnight with chicken anti-GFP primary antibody (Abcam, ab13970; 1:1000) diluted in PBST. Then, we washed the samples with PBST (3 times, 40 minutes each) and incubated them overnight with secondary antibodies diluted in PBST. Finally, we washed the samples with PBST (3 times, 40 minutes each), placed them on glass slides, and mounted coverslips with ProLong Gold (Invitrogen, 36934).

Analysis of ear skin whole mounts

To assess branch points per millimeter of vessel, we measured the length of the lymphatic network and then manually counted the number of lymphatic branch points. To measure vessel diameters, we placed a 15 × 15 grid over images and measured the diameter of vessels located on intersecting grid lines. To assess lymphatic valves per millimeter of vessel, we measured the length of the lymphatic network and then manually counted the number of lymphatic valves. We used NIS Elements software (v5.30.02) to manually count branch points and valves and to measure vessel lengths and diameters.

Lung inflation apparatus assembly

A 20 mL Luer-Lok tip syringe was attached to a ring stand 25 cm above the animal dissection platform. A stopcock was twisted into the bottom of the Leuer-Lok tip syringe, and rubber tubing was connected to the stopcock. A blunt-ended needle was coupled to an adapter joined to the rubber tubing.

Lung inflation, sectioning, staining, and analysis

Mice were euthanized and then perfused with heparin/DPBS pH 7.2 to remove red blood cells. The blunt-ended needle joined to the lung inflation apparatus was inserted into a small incision in the trachea and secured with sutures. The stopcock was opened, allowing 4% PFA to inflate the lungs for 5 minutes. The lungs were then collected, and the suture was tied at the end of the trachea. The lungs were placed in a 50 mL conical tube containing 4% PFA, and the sutures were held in place by the threads of the cap and tube. The 50 mL conical tube was then inverted, allowing the lungs to be fully submerged while fixing overnight at 4°C. Samples were washed with 1× PBS 3 times (5 minutes each) and placed in cryoprotect solution (20% sucrose + 2% polyvinylpyrrolidone in 1× PBS) overnight at 4°C. The next day, samples were transferred to a 50 mL conical tube containing prewarmed (54°C) gelatin embedding medium (8% powdered gelatin + 20% sucrose + 2% polyvinylpyrrolidone in 1× PBS) and stored overnight at 54°C. Samples were then placed in cryomolds containing gelatin embedding medium, cooled overnight in the refrigerator, and stored in a –80°C freezer. We used a cryostat to cut 100-μm-thick sections, which we placed on glass slides. The slides were then placed in a 37°C incubator for 10 minutes. The tissues were transferred to a 24-well plate containing water. The samples were washed 3 times (5 minutes each) and then permeabilized with PBST for 1 hour at room temperature. We then followed our protocol for whole-mount immunofluorescent staining of mouse ear skin. We used a goat anti-Vegfr3 antibody (R&D Systems, AF745; 1:250). Samples were transferred to glass slides and cover slips mounted with ProLong Gold (Invitrogen, 36934). Images were captured with a 4× objective, and we used NIS Elements software (v5.30.02) to measure vessel diameters. Briefly, a 50 μm × 50 μm grid was placed over the images, and vessel diameters were measured at sites where there were intersecting gridlines.

Analysis of tissue sections

Ear skin was fixed overnight in 4% PFA, washed with PBS (3 times, 5 minutes each), and submitted to the histology core for processing. A microtome was used to cut 5-μm-thick sections, which were then placed on microscope slides. Slides were deparaffinized with xylene and rehydrated through a descending ethanol series. Antigen retrieval was performed by heating slides in 0.01 M citric acid (pH 6.0) in a pressure cooker. Nonspecific binding was blocked by incubating slides with Tris-buffered saline with 0.2% Tween 20 (TBST) plus 3% donkey serum for 30 minutes. Slides were then incubated overnight with goat anti-Prox1 (R&D Systems, AF2727; 1:250) and rabbit anti-Lyve1 (Abcam, ab14917; 1:250) primary antibodies diluted in TBST plus 5% BSA. Slides were washed with TBST and then incubated with fluorophore-conjugated secondary antibodies diluted in TBST plus 5% BSA. Following washes with TBST, coverslips were mounted with ProLong Gold plus DAPI (Invitrogen, 36935). Four to 12 lymphatics were analyzed per sample. We used NIS Elements software (v5.30.02) to measure the length of Lyve1-positive lymphatic vessels and count the number of Prox1-positive nuclei per lymphatic vessel.

Intralobular lymphatic vessel imaging

Flt4CreERT2;Kras+/G12D mice were crossed with the Rosa26mT/mG reporter line to visualize all lymphatic vessels. Tamoxifen (100 μg) was administered subcutaneously to newborn pups on P0 and P2 to activate Cre recombinase in VEGFR3-expressing cells and tissues were analyzed 3 weeks later. Mice were euthanized on P21 and whole lungs were excised from the chest cavity with the trachea and heart attached. The lungs were then placed in a dish containing a silicon bottom (made of Sylgard 170) and minuten pins were used to position the tissue with the ventral side facing up while bathed in ice-cold PBS. Images were immediately taken of the collecting lymphatic vessels entering each lobe on a Zeiss AxioZoom V16 microscope.

Statistics

Data were analyzed using GraphPad Prism statistical analysis software (v9.5.1). A Shapiro-Wilk normality test was performed to assess whether the data were normally distributed. If the data were normally distributed, 2-tailed, unpaired Student’s t tests were performed to test means for significance. When the data were not normally distributed, the Mann-Whitney test (2-tailed) was conducted. All results are expressed as mean ± SEM. The number of mice in each group is indicated in the figure legends (n = number of mice). A P value of less than 0.05 was considered significant.

Study approval

The animal experiments described in this study were carried out in accordance with an animal protocol approved by the Institutional Animal Care and Use Committee of UT Southwestern Medical Center.

Data availability

All sequencing data have been deposited in the NCBI Gene Expression Omnibus (GEO) under accession number GSE272843 and are publicly available as of the date of publication. Raw data can be found in the Supporting Data Values file.