Methods

Wildtype C57BL/6J mice were purchased from Janvier Labs. The transgenic lines B6 _Evi/Wls_fl/fl × Cdh5-CreERT2 and B6 _Rspo3_fl/fl × Cdh5-CreERT2 were described previously82,83. The double mutant line was generated by crossing _Evi/Wls_fl/fl and _Rspo3_fl/fl with B6 Cdh5-CreERT2 line. Cre-positive mice used in this study were homozygous for the floxed alleles and heterozygous for Cre. Cre-negative littermates that were homozygous for the floxed alleles served as controls. All animals were housed and maintained in pathogen-free animal facilities of the German Cancer Research Centre. In brief, all mice were housed on a temperature controlled 12 h light/12 h dark cycle and had free access to water and food ad libitum. In vivo experiments were performed under the guidelines of the institutional and governmental regulation approved by the regional council Karlsruhe (DKFZ305, DKFZ370, DKFZ434, G-255/18, G-107/18, G-251/20 and G-36/22).

For the samples of portosystemic shunt, liver sections were shared and collected from Heikenwalder lab29 and Kiefer lab. For the murine liver samples directly exposed to normobaric hypoxia (10% O2), the samples were obtained from Korff lab40. In brief, mice were exposed to normobaric hypoxia for 3 or 7 days in a hypoxia chamber (Biospherix) with free access to water and food ad libitum with 10% oxygen and 90% nitrogen.

Patient specimen

Patient samples were obtained under licence 2022-140_1_S_KH approved by the local ethics committee of TUM/MRI under the legal guidelines of the Government of Upper Bavaria (BayKrG Art. 27 Abs. 4), with a waiver of informed consent due to retrospective and anonymised secondary use of clinical specimens. The samples were subjected to immunostaining for GS.

Tamoxifen treatment

To induce genetic recombination, 6–10 weeks old male and female mice were intraperitoneally (i.p.) injected with 2 mg/kg tamoxifen dissolved in 50 μl of corn oil for 5 consecutive days followed by a two-week washout period. For long-term tracing experiments, mice were sacrificed 3, 6 and 12 months after the initial genetic recombination. To induce genetic recombination in vitro, hydroxytamoxifen (4-OHT; Sigma) was dissolved in DMSO (Applichem) at a stock concentration of 20 mM and primary LSEC were treated with a final concentration of 2 μM of 4-OHT for two days with daily medium changes. Both female and male mice were included in our in vitro study.

Partial hepatectomy

Two-third PHx was performed as described previously22. In brief, mice were anaesthetised with a mixture of ketamine (Bela-Pharm, 100 mg/kg body weight) and xylazine (Bayer, 10 mg/kg body weight) by i.p. injection. Next, mice were put on a temperature-controlled pad and eyes were protected from drying by using an eye cream (Bayer). Mice were laid on their back and skin and peritoneum were incised to allow full access to the liver. Subsequently, the left lateral lobe and the median lobe of the liver were ligated using 4-0 silk sutures (Ethicon) and resected. For one-third PHx, the left lateral lobe of the liver was resected. Metamizole (Serumwerk Bernburg) was used as a post-surgical analgesic treatment (200 mg/kg body weight for subcutaneous injection and 800 mg/kg body weight for ad libitum) during the first 48 h post-surgery. Mice were euthanised at the indicated time points using cervical dislocation. Both male and female mice were incorporated in the study.

Ex vivo perfusion

Mice were first sacrificed by cervical dislocation and livers were exposed for ex vivo perfusion. The portal vein was cannulated with 27 G winged infusion system (B.Braun) and connected to a peristaltic pump (Ismatec). Subsequently, livers were perfused using 37 °C preheated Krebs–Henseleit saline solution (HiMedia) supplemented with sodium bicarbonate (Sigma) for 30 min. To mimic the physiological flow, the flow rate of the pump was set to 4 ml/min. In parallel, 8 ml/min of flow rate was applied to mimic the accelerated flow23,38. Both male and female mice were included in the study.

Pharmacological in vivo treatments

To induce vasodilation, female and male mice were treated with Riociguat (AdooQ Bioscience). The drug was dissolved in DMA (Sigma) in a concentration of 100 mg/ml and stored at −20 °C. Before the experiment, Riociguat/DMA solution was further diluted in 20% 2-Hydroxypropyl-β-cyclodextrin/H2O solution (Sigma) to reach a final concentration of 2.5 mg/ml37. Mice were treated with Riociguat (10 mg/kg body weight) or vehicle control twice a day via oral gavage for two days.

To induce hypoxia in vivo, female and male mice were treated with Roxadustat (Selleckchem), a HIF prolyl hydroxylase Inhibitor. Roxadustat was dissolved in DMSO in a concentration of 50 mg/ml and stored at −20 °C. The working solution was freshly prepared by further diluting the solution with a sterile PBS (Sigma) in a concentration of 1 mg/ml84. For the pulse treatment, wild-type mice were injected i.p. with a dose of 10 mg/kg and sacrificed 4 h after the injection. For the long-term treatment, wild-type mice received daily i.p. injections for five consecutive days at a dose 10 mg/kg body weight.

To reduce hepatic oxidative stress, NAC (Sigma) was added to the drinking water of wild-type animals at a concentration of 1 g/L for two weeks85. Both male and female mice were incorporated in the study.

Isolation of hepatocytes and liver sinusoidal endothelial cells

The optimised protocol of LSEC isolation with a minimum shear activation was adapted from the previous studies22. In brief, both female and male mice were first sacrificed by cervical dislocation and livers were exposed for in situ perfusion. For this, the vena cava was fixed using a 27 G winged infusion set connected to the tubing of a peristaltic pump. The liver was perfused with 37 °C pre-warmed liver perfusion medium (Gibco) at 4 ml/min for 2 min, followed by 37 °C pre-warmed liver digestion medium (Gibco) supplemented with 20 μg/ml Liberase TM (Roche) at 3 ml/min for 5 min. Perfused livers were explanted and placed into a Petri dish with a pre-warmed RPMI medium (Gibco). After removing the liver capsule membrane, the tissue was dissociated by gentle shaking in a final volume of 40 ml of RPMI.

For hepatocyte collection, the single cell suspension was filtered through a 100μm cell strainer (BD Bioscience) and centrifuged once at 50 × g and 4 °C for 5 min. The resulting cell pellet contained hepatocytes, which were used for further experiments.

For the isolation of non-parenchymal hepatic cells (NPC), the supernatant was collected and centrifuged for an additional 5 min at 50 × g and 4 °C to remove excessive parenchymal cells. The supernatant was collected, and the final NPC population was obtained by centrifugation at 300 × g and 4 °C for 10 min. Red blood cells were lysed by incubating the isolated cells in ACK buffer for 5 min. To purify LSEC, total NPC were resuspended in 200 μl of FACS buffer containing 5% FCS and incubated with 20 μl of mouse CD146 MicroBeads (Miltenyi Biotec) for 20 min on ice. The cells were washed once, resuspended in 1 ml of FACS buffer and loaded onto a LS column (Miltenyi Biotec) and cells were purified according to the manufacturer’s instructions. The isolated LSEC were then processed for flow cytometry, frozen for further molecular analyses or taken into cell culture as freshly isolated primary LSEC.

Flow cytometry and spatial cell sorting

For hepatocyte cell sorting, enriched isolated hepatocytes (see detailed protocol above) were resuspended in FACS buffer and sorted using BD FACS Aria I and II (BD Bioscience). Hepatocytes were gated based on (1) forward scatter area (FSC-A) and side scatter area (SSC-A) to exclude NPC and cell debris, (2) forward scatter width and forward scatter height (FSC-H), as well as side scatter width and side scatter height to exclude cell doublets, and (3) FxCycle Violet staining at a dilution of 1:1000 (Invitrogen) to exclude dead or damaged cells. For spatial sorting, hepatocytes were stained with CD73-FITC (BioLegend) and CD324-PE (BioLegend) at a dilution of 1:200 on ice for 20 min. Cells were washed twice and resuspended in 1 ml FACS buffer. To obtain hepatocytes from different zones, the previously described FACS gating strategy was applied (Supplementary Fig. 6a, m)53.

For FACS sorting of LSEC, the NPC cell population was stained with CD45-FITC (BD Bioscience), CD31-APC or PE (BD Bioscience), and CD146-PE-Cy7 (Biolegend) at a dilution of 1:200 for 20 min. Stained cells were washed and resuspended in 1 ml FACS buffer, and cell populations were sorted using a BD FACSAria cell sorter platform. LSEC were isolated based on (1) FSC-A and SSC-A to exclude cell debris and cell clumps, (2) FSC-A and FSC-H to exclude cell doublets, (3) FxCycle Violet stain at a dilution of 1:1000 to exclude dead or damaged cells, and (4) CD45 negative, and CD31 and CD146 double positive staining to exclude contaminating cell types. To perform spatial sorting of LSEC, a previously described strategy was adapted22. In brief, the NPC cell population was additionally stained with CD117 (c-Kit)-APC (BioLegend) at a dilution of 1:200, and LSEC from the distinct zones were sorted based on their c-Kit expression profile (Supplementary Fig. 6b, n).

For FACS-based cell size analysis, spatial sorting of LSEC was performed to distinguish portal, midzonal, and central LSEC. FSC-A was further recorded according to the previously described FACS gating strategy (Supplementary Fig. 9j). The geometric mean of FSC-A was calculated, and the abundance of LSEC according to its FSC-A value was plotted on FlowJo (BD Bioscience, v.10).

Cell culture and in vitro studies

Freshly isolated primary LSEC (see above for detailed protocol) were seeded into 6-well plates at a density of 1.0 × 106 cells per well. The cells were cultured in a complete Endothelial cell medium (Cell biologics) with full supplementation according to the manufacturer’s instructions and maintained at 37 °C under 5% CO2.

Manipulation of hypoxia in vitro

Freshly isolated primary LSEC obtained from wild-type mice were cultured in a 1% oxygen cell culture chamber at 37 °C for 24 h and cells were collected for further analysis. To pharmacologically mimic hypoxic conditions, cells were treated with 100 and 200 μM of cobalt-(II)-hexahydrate (CoCl2; Sigma) for 24 h. Cells were then washed and processed for further analysis.

Manipulation of shear stress in vitro

To induce endothelial shear response in vitro86, freshly isolated wild-type primary LSEC were placed on a temperature-controlled orbital shaker (Edmund Buehler) set to 37 °C and were shaken at 90 rpm mimicking low shear conditions with a shear of 7 dynes/cm2 or at 200 rpm simulating high shear conditions with a shear of 14 dynes/cm2. These cells were shaken for 12 h and samples were collected for secondary analysis.

To discriminate between laminar shear and oscillatory shear, shear stress was induced using the in vitro fertilization (IVF) dish (Falcon) under an orbital shaker. Shear force at the periphery area of IVF dish mimics laminar shear. Whereas the shear force at the centre area of the IVF dish mimics oscillatory shear. Eight wild-type mice were pooled to isolate 40 million primary LSEC. Next, two million LSEC were seeded either on the central area or the peripheral area of the IVF dish overnight. For the dose-effect experiment, a different magnitude of shear stress was applied using an orbital shaker (0, 7, 14 dynes/cm2) for 12 h. To capture the temporal kinetics of the LSEC transcriptome, two million LSEC (pooled from eight mice) were seeded on the central area of the IVF dish overnight and oscillating shear stress was applied using an orbital shaker at 14 dynes/cm2 for different time points (30 min, 1, 6 and 12 h).

To further evaluate the shear response of primary LSEC, a conventional IBIDI microfluidic pump system (IBIDI) was used. In brief, 1.0 × 106 primary LSEC were seeded on 0.1% Gelatin-coated IBIDI Luer 0.4 mm slide (IBIDI) overnight. The slides were then subjected to either static control, laminar shear, or oscillatory shear. To induce laminar shear, a single IBIDI pump unit was acquired and set to provide a shear of approximately 14 dynes/cm² overnight. Conversely, two IBIDI pump units were used to induce oscillatory shear stress with a magnitude of approximately 14 dynes/cm² overnight. After the experiments, cell morphology was confirmed on a microscope and cells were directly lysed for mRNA extraction.

Pharmacological manipulation of Wnt signalling pathway in vitro

To activate canonical Wnt signalling, freshly isolated primary LSEC were treated with 3 μM CHIR99021 (Selleckchem), a GSK3β inhibitor, for 24 h. Conversely, the canonical Wnt signalling pathway was inhibited using a porcupine inhibitor LGK974 (Selleckchem) at a concentration of 10 nM for 24 h. To induce the non-canonical Wnt signalling pathway, the cells were treated with recombinant WNT5A (R&D systems) at 10, 100 or 200 ng/ml for 24 h.

To mimic the central venous milieu, freshly isolated primary LSEC were treated for 24 h with 200 ng/ml of recombinant WNT2 (Cusabio) or 100 ng/ml of recombinant WNT9B (R&D).

TCF reporter assay

To test for the secretion of bioactive Wnt, wild-type primary LSEC were cultured in serum free EC medium and subjected to 12 h of shear activation. The conditioned medium (CM) was collected. 1.0 × 105 HEK293T cells transfected with the TCF4/Wnt-mCherry construct were seeded in a six-well dish and cultured with complete DMEM (i.e., supplemented with 10% FCS and 1% penicillin/streptomycin). After one day, the complete DMEM medium was replaced with a 1:1 mixture of complete DMEM medium and primary LSEC CM from the different shear conditions. This mixture was further supplemented with either 100 ng/ml of canonical WNT3A (PeProtech) or 1 μM of CHIR99021 (Selleckchem) to serve as positive controls. HEK293T reporter cells were incubated for 12 h and the TCF4/Wnt activity was recorded using flow cytometry. To quantitatively assess the Wnt-activation, the geometric mean intensity of TCF4-induced mCherry was calculated using FlowJo (BD Bioscience, v.10).

Live cell calcium imaging

Primary LSEC were isolated from the double mutant and control mice. 1.0 × 106 cells were plated in a 6-well plate and genetic recombination was induced as mentioned above. For flow-dependent study, primary LSEC were seeded on the periphery area of the IVF and laminar shear was applied using an orbital shaker. Genetically modified primary LSEC and control LSEC were incubated with 3 μM of Cal-520 (AAT Bioquest) in EC medium in a humidified incubator at 37 °C for 30 min. For a control experiment, cells were co-incubated with 100 μM of Carbenoxolone (CBX, Sigma-Aldrich) for 30 min. The cells were washed using Hanks Buffer (Gibco) and submerged in calcium-depleted Hanks Buffer (Gibco). Next, the plate was placed on Zeiss Widefield Cell Observer. Internal calcium storage was depleted using 1 μM of Thapsigargin (Sigma) and the image was acquired for the baseline control. Subsequently, 1 mM of CaCl2 was added to the cells submerged with the calcium-depleted Hanks Buffer and the image was acquired every second for 5 min. Next, 1 μM of ionomycin (Sigma-Aldrich) was added to the solution as a positive control.

RNA extraction, cDNA synthesis, and real-time quantitative RT-PCR

For RNA extraction of whole liver lysates, liver pieces were homogenised using TissueLyzer (Retsch) in TRIzol (Sigma) and RNA was purified using chloroform (Carl Roth) and isopropanol (Sigma). RNA extraction from cell pellets (sorted or cultured cells) was performed using the GeneElute mammalian Total RNA Isolation Kit (Sigma) according to the manufacturer’s instructions. 500 to 1000 ng of total mRNA was reverse transcribed into cDNA using the QuantiTect Reverse transcription kit (Qiagen) following the manufacturer’s instructions. Gene expression analysis was performed by quantitative PCR using TaqMan reactions (Thermo Fisher Scientific) and StepOnePlus (Applied Biosystems). Gene expression levels were assessed using the Ct-method and normalised to the expression of the respective control gene, resulting in ΔCt-values. Relative gene expression was assessed by normalising ΔCt-values of individual samples to the average control ΔCt-value, resulting in ΔΔCt-values. Relative fold changes to control were then calculated as 2−ΔΔCt. The following taqman probes (Thermo Fisher Scientific) were used: Klf2, Mm01244979_g1; Angpt2, Mm00545822_m1; Wnt2, Mm00470018_m1; Wnt9b, Mm00457102_m1; Gja1, Mm00439105_m1; Gja4, Mm00433610_s1; Klf4, Mm00516104_m1; Nos3, Mm00435217_m1; Pecam1, Mm01242584_m1; Actb, Mm00607939_s1; Sox9, Mm00448840_m1; Mfsd2a, Mm01192208_m1; Igfbp2, Mm00492632_m1; Hamp2, Mm00842044_g1; Axin2, Mm00443610_m1; Lgr5, Mm00438890_m1; Fzd4, Mm00433382_m1; Lrp6, Mm00999795_m1; Rspo3, Mm01188251_m1; Wls, Mm00509695_m1; Hif1a, Mm00468869_m1; Epas1, Mm01236112_m1; Slc2a1, Mm00441480_m1.

Protein extraction and immunoblotting

Proteins were isolated using RIPA buffer (Thermo Scientific) for cells and Trizol for liver tissues. Proteins were quantified using BCA (Thermo Scientific) and adjusted to the same concentration. Next, the samples were heated at 95 °C for 5 min, and subjected to 10% SDS-polyacrylamide gel electrophoresis. Proteins were transferred to PVDF membranes and were further blocked with 3% BSA (Gerbu) in TBST buffer. Following primary antibodies were used for immunoblotting: rabbit anti-mouse GJA1 (1:2000, Abcam), rabbit anti-mouse GJA4 (1:1000, Invitrogen), rat anti-mouse activated-CD29 (1:1000, BD), goat anti-mouse VEGFR3 (1:1000, R&D), rabbit anti-mouse Phospho-VEGFR3 Tyr1230 (1:1000, Affinity Biosciences), rabbit anti-mouse KLF2 (1:1000, Milipore), goat anti-mouse KLF4 (1:1000, R&D), mouse anti-mouse HIF1A (1:1000, Novus Biologicals), and rabbit anti-mouse Tubulin (1:1000, Cell Signaling). Consequently, membranes were washed with TBST, incubated with HRP-conjugated secondary antibodies: goat anti-rabbit HRP (1:5000, DAKO) and rabbit anti-goat HRP (1:5000, DAKO), developed with Western Blotting Lumino Reagent (Santa Cruz), and images acquired by Amersham Imager 600 (GE Healthcare). For reblotting, membranes were stripped using Re-blot strong solution (Millipore). Quantification of immunoblotting was performed using FIJI software (ImageJ, 1.54 g).

Sample preparation of metabolomic screening

Blood plasma and liver pieces (lateral lobe) from 3-month-old mutant and control female and male mice were collected. The frozen liver pieces were snap frozen in liquid nitrogen and pulverised using a cooled ball mill (components were pre-cooled in liquid nitrogen). For the semi-targeted metabolite screening, 50 mg frozen ground liver tissue or 50 µL plasma were processed in 360 µL of 100% MeOH for 15 min at 70 °C with vigorous shaking. As an internal standard, 20 µL of Ribitol (Sigma) at a concentration of 0.2 mg/ml was added to each sample. The samples were then incubated in 200 µL of chloroform and shaken for 5 min at 37 °C. To separate polar and organic phases, 400 µL water was added and samples were centrifuged for 10 min at 11,000 × g. For the derivatisation, 700 µL of the polar (upper) phase was transferred to a fresh tube and dried in a speed-vac (vacuum concentrator) without heating.

Derivatisation (methoximation and silylation)

Pellets of the aqueous phase after extraction were redissolved in 20 µL of methoximation reagent containing 20 mg/ml methoxyamine hydrochloride (Sigma) in pyridine (Sigma) and incubated for 2 h at 37 °C with shaking. For silylation, 32.2 µL of N-Methyl-N-(trimethylsilyl)trifluoroacetamide (Sigma) and 2.8 µL of Alkane Standard Mixture (50 mg/ml C10–C40; Fluka) were added to each sample. After 30 min of incubation at 37 °C, the samples were transferred to glass vials for GC/MS analysis.

Gas Chromatography/Mass Spectrometry (GC/MS) analysis

GC/MS was performed using GC/MS-QP2010 Plus (Shimadzu) fitted with a Zebron ZB 5MS column (Phenomenex; 30 m × 0.25 mm × 0.25 µm). The GC was operated with an injection temperature of 230 °C and 1 µL sample was injected with split mode (1:5). The GC temperature program started with a 1 min hold at 40 °C followed by a 6 °C/min ramp to 210 °C, a 20 °C/min ramp to 330 °C and a bake-out for 5 min at 330 °C using Helium as the carrier gas with constant linear velocity. The MS was operated with ion source and interface temperatures of 250 °C, a solvent cut time of 7 min and a scan range (m/z) of 40–700 with an event time of 0.2 s. The GC-MS solution software (Shimadzu) was used for data processing.

For the measurement of cations, 30 µl extract (aqueous extract prepared for GC/MS analysis) was diluted with 270 µl ultra-pure water. Cations were separated isocratically with 30 mM methansulfonic acid for 27 min on an IonPac CS16 column (2 mm, Thermo Scientific) connected to an ICS-1000 system (Dionex) at a flow rate of 0.36 ml/min and quantified by conductivity detection after anion suppression (CERS-500 2 mm, suppressor current 43 mA) using ultra-pure standards (Sigma). Column temperature was set to 43 °C. Metabolomic data can be accessed on Supplementary Data 4.

Sample preparation for proteomic analysis

Cell pellets were dissolved in an appropriate volume of freshly prepared 0.1% Rapigest (100 mM ammonium bicarbonate). The samples were sonicated at 10% amplitude on 50% Duty cycle with 20 cycles for 80 s. The cell lysates were centrifuged at 15,000 × g at 4 °C for 30 min and the protein concentration of the supernatants was measured using BCA protein assay kit (Pierce). A sufficient protein amount from each sample was reduced by incubating with 10 mM DTT at 56 °C for 30 min. The samples were then alkylated using chloroacetamide to a final concentration of 50 mM and incubated for 30 min at RT. The alkylation was quenched by adding DTT to a final concentration of 10 mM. Proteins were digested using trypsin–LysC at a 1:50 ratio and incubated overnight at 4 °C. The samples were acidified by adding an appropriate amount of 10% TFA and the acidity (pH <2) of the solution was checked. The samples were then incubated at 37 °C for 30 min and centrifuged for 30 min at maximum speed. The supernatants were transferred to fresh tubes and stored at −20 °C until the time of the analysis. Both female and male mice were incorporated in the proteomic screening.

Liquid chromatography with tandem mass spectrometry (LC MS/MS) analysis

Peptides were separated using an EASY nanoLC 1200 System (Thermo Scientific), equipped with Acclaim PepMap 100, C18 trapping column (Thermo Scientific, 100 μm × 2 cm, 5 μm, 100 Å) and Acclaim PepMap RSCL 2 mM C18 analytical column (Thermo Scientific, 75 mm × 50 cm). The analytical column was heated to 45 °C in a MonoSLEEEVE column oven (Analytical Sales and Services). Sample loading on the trap column was performed under a constant pressure of 800 bar for a total volume of 22 μl of solvent A (0.1% Formic acid). Peptides were eluted from the trapping column with the following gradient at a flow rate of 300 nl/min: 0 to 4 min linear gradient from 3 to 8% solvent B (80% Acetonitrile and 0.1% Formic acid), 4 to 6 min linear gradient 8 to 10% B, 6 to 74 min linear gradient from 10 to 32% B, 74–86 min linear gradient from 32 to 50% B, 86 to 87 min gradient increase from 50 to 100% B, 87 to 94 min the gradient was kept to 100% B, 94 to 95 min linear gradient from 100 to 3% B, 95 to 105 min gradient was kept at 3% B.

The analytical column was coupled to a Silica PicoTip Emitter (Scientific Instrument Services) via a zero dead volume connection. Peptides were ionised using a Nanospray Flex nano spray source (Thermo Scientific) at a 2 kV voltage in the QE–HF mass spectrometer (Thermo Scientific). The ion transfer tube temperature was set to 275 °C. The following settings were used for data acquisition: Full scan MS1 acquisition was acquired within 350–1500 m/z scan range and 120,000 resolution. The maximum injection time (IT) was set to 32 ms and automatic gain control (AGC) to 3 × 106 ions. Spectrum data type was set to profile. The top 20 most abundant precursor ions were selected for fragmentation. Ions with unassigned charges and charges of 1 or >5 were excluded. Dynamic exclusion was set to 40 s with a mass tolerance of ±10 ppm. For the MS2 scans, quadrupole was used with an isolation window of 2.0 m/z. For peptide fragmentation, higher-energy collisional dissociation was used at 26%. MS2 scans were acquired with an AGC target of 1 × 105 or maximum IT of 50 ms and scan range window between 200 and 2000 m/z. MS2 spectra were acquired in centroid mode.

Analysis of raw MS data

The raw data were searched using Maxquant (v 1.6.10.43) against the mouse canonical and reviewed database (Mus musculus, April, 2019). Enzyme digestion was set to trypsin allowing for a maximum of up to 2 missed-cleavages. Oxidation (M), Acetyl (Protein N-term), and deamidation (NQ) were set as variable modifications and carbamidomethyl (C) as fixed. Minimum unique peptides for protein identification were set to 1, while peptide and protein hits were filtered at a false discovery rate (FDR) of 1% with a minimal peptide length of 7 amino acids. The reversed sequences of the mouse canonical database were used as a decoy database for FDR calculation. Second peptide search option was enabled. Match between runs option was enabled with a matching time window of 0.4 min and the alignment time window set to 20 min. Dependent peptides option was deactivated. iBAQ calculation was also enabled. Label free quantification option was also enabled and LFQ min.ratio.count option was changed to 1. All other Maxquant options were left to their default settings. LFQ-intensity of hepatocytes and LSEC are available on Supplementary Data 4.

Further gene ontology enrichment and pathway analysis on MS data were performed as previously described using LFQ-Anaylst (v 1.2.4)87. In brief, the pathway analysis was computed based on the LFQ counts. Differentially expressed proteins were calculated using limma from R for the pair-wise comparison (adjusted p < 0.2 with Benjamini–Hochberg method & log2 fold change > 1.0), consequently, pathway analysis was performed.

Immunofluorescence imaging

Liver lobes were fixed in 4% paraformaldehyde solution overnight and dehydrated in 30% sucrose prior to embedding in the OCT compound (Sakura Finetek). Cryosections (10 μm) were cut using a CM3050S cryostat (Leica) and attached to Superfrost Plus slides (VWR). Antigen retrieval was performed for APC, FZD4, and LRP6 staining by incubating the tissue section with Proteinase K (20 µg/ml, Gerbu Biotechnik) for 15 min at 37 °C. The sections were then submerged in a blocking buffer composed of PBS with 0.3% TritonX-100 (Sigma-Aldrich) and 10% FCS (Gibco) for 1 h in RT. Subsequently, the primary antibodies were diluted in the blocking solution and incubated overnight at 4 °C. Following primary antibodies were used for staining: rabbit anti-mouse GS (1:1000, Abcam); rabbit anti-mouse Cytochrome P450 2E1 (CYP2E1, 1:500, Abcam); goat anti-mouse VE-Cadherin (1:200, R&D); goat anti-mouse ICAM1 (1:200, R&D); mouse anti-mouse eNOS (1:200, Cell Signaling); goat anti-mouse Collagen IV (COLIV, 1:200, R&D Systems); rabbit anti-mouse GJA4 (1:200, Invitrogen); rabbit anti-mouse APC (1:100, Novus); goat anti-mouse FZD4 (1:100, Novus); goat anti-mouse LRP6 (1:100, Novus); rat anti-mouse KI67 (1:200, Invitrogen); rabbit anti-mouse HNF4A (1:200, Abcam); rabbit anti-mouse GJA4 (1:1000, Invitrogen). The sections were then washed three times using a wash buffer composed of PBS with 0.2% TritonX-100 and 5% FCS. Subsequently, sections were incubated with secondary antibody solution diluted in washing buffer at a dilution of 1:400 for 2 h at RT. The following secondary antibodies were used for staining: Alexa Flour 546 Phalloidin, Alexa Fluor 647 donkey anti-rabbit IgG (H + L); Alexa Fluor 568 donkey anti-goat IgG (H + L) (Thermo Fisher Scientific); Alexa Fluor 488 donkey anti-rat IgG (H + L). Stained sections were counterstained with 1:2000 Hoechst 33342 (Sigma-Aldrich) and mounted with Fluorescence Mounting Medium (Agilent Dako). Images were acquired with either Leica TCS SP5II (Leica) or Leica TCS SP8 (Leica) confocal microscope and image analysis was performed using Fiji software (ImageJ, 1.53t). For high resolution imaging, images were acquired using LSM 980 Airyscan (Zeiss) with Airyscan 2 modes focussing on resolution and sensitivity. The images were processed using ZEN blue (Zeiss).

For staining of hypoxic areas, pimonidazole staining (Hypoxyprobe) was applied according to the manufacturer’s protocol. In brief, mice were intravenously injected with pimonidazole (60 mg/kg) and organs were harvested after 1 h of exposure. Tissue sections were prepared as described above and immunostaining was carried out according to the manufacturer’s instructions.

Image analysis was performed using Fiji software (ImageJ, 1.53t). After region-of-interest (ROI) selection, GS, COLIV, LRP6, APC, FZD4, KI67, HNF4A, DAPI channels were masked using thresholding. For liver zonation, the percentage of GS+ area within the ROI was calculated. For vascular Wnt zonation, overlapping of FZD4+ and COLIV+, or LRP6+ and COLIV+, or APC+ and COLIV+ LSEC within the ROI were quantified. To quantify regenerating hepatocytes, KI67+/HNF4A+/DAPI+ triple-positive hepatocytes were considered proliferating hepatocytes. For analysing zone-specific hepatocyte proliferation, ROIs reflecting portal, midzonal, and central regions were defined based on the spatial coordination of GS+ hepatocytes and bile ducts. Subsequently, KI67+/HNF4A+/DAPI+ triple-positive hepatocytes in each of the three zones were quantified.

Proximity ligation assay

The proximity ligation assays (PLA) were performed using the Duolink Proximity Ligation Kit (Sigma) and the manufacturer’s protocol was followed with minor modifications. Cryosections were prepared (7 μm) from the PFA-fixed liver samples of both female and male mice and blocked using the Duolink Blocking Solution at 37 °C for 60 min in a humid environment. Primary antibodies were diluted in Duolink Antibody Diluent in a 1:100 ratio and incubated at 37 °C overnight in a humid environment. The following antibodies were used for PLA: Rat anti-mouse ITGB1 (1:100, BD), Goat anti-mouse VEGFR3 (1:100, R&D), Rabbit anti-mouse Phosphotyrosine (1:100, Abnova), and Rabbit anti-mouse Phosphoserine (1:100, Abnova). The antibodies were further labelled with Duolink anti-rabbit MINUS probe (Sigma) and Duolink anti-goat PLUS probe (Sigma), while ITGB1 antibody was conjugated with Duolink in situ Probemaker PLUS (Sigma). PLA probe incubation and ligation were performed according to the manufacturer’s protocol. For signal amplification, Duolink In Situ Detection Reagents Orange (Sigma) was applied. For counterstaining, goat anti-mouse VE-Cadherin (1:200, R&D) or rat anti-mouse ITGB1 (1:200, BD) was used. Furthermore, rabbit anti-Glutamine synthetase (1:1000, Abcam) was pre-labelled with Alexa Fluor 647 using Zenon IgG labelling reagent (Thermo Scientific) and stained to distinguish liver zonation.

In situ hybridisation assay

The in situ hybridisation assays were performed using the ViewRNA ISH Tissue Assay Core Kit (Thermo Fisher Scientific). The manufacturer’s protocol was employed with minor modifications. In brief, cryosections (7 μm) were fixed with PFA overnight at 4 °C followed by dehydration in ethanol, baked for 1 h at 60 °C, boiled for 20 min in pre-treatment solution, and further digested for 20 min in protease solution at 42 °C. Subsequently, the hybridisation was performed with Gja1 (VB6-17027-VT) probe (Thermo Fisher Scientific) for 2 h at 42 °C. The signals were then pre-amplified and amplified according to the protocol at 42 °C. Next, Gja1-labelled sections were conjugated to the fast red substrate (Thermo Fisher Scientific) for visualisation. Consequently, goat anti-mouse VE-Cadherin (1:200, R&D) and rabbit anti-Glutamine synthetase (1:1000, Abcam) were treated for counterstaining.

Shear stress–induced cell alignment

Phalloidin staining was performed on either static or shear activated primary LSEC isolated from mutant and control mice. First, 1.0 × 106 of primary LSEC were seeded and were grown until they formed a monolayer. Next, in vitro gene deletion was induced using 4-OHT as described above. Cells were then placed in an orbital shaker for shear activation at 14 dynes/cm2 for 12 h. For peristaltic pump-induced shear stress, primary LSEC were seeded on 0.1% Gelatin-coated IBIDI Luer 0.4 mm slide (IBIDI) overnight. Subsequently, in vitro gene deletion was performed for 48 h. To avoid cell detachment, a transient pulse of low shear stress (3 dynes/cm²) was applied for 30 minutes, and ultimately, the LSEC were exposed to high laminar shear stress (11 dynes/cm²) for 24 h. Next, cells were fixed in 4% PFA for 30 min and blocked using a blocking solution composed of 0.3% TritonX-100 and 10% FCS mixed in PBS. Subsequently, phalloidin (Invitrogen) with cell nucleus staining was performed for 1 h in RT and cells were washed and submerged in PBS. Imaging was performed on Zeiss Widefield Axio CellObserver (Zeiss). Six random peripheral areas of each well were selected and imaged. The images were then quantified and analysed using Fiji (ImageJ, 1.53t). In brief, the stress fibres marked by F-actin were segmented using OrientationJ (ImageJ, 1.53t). Then the direction and angle of segmented filaments were calculated through Directionality (ImageJ, 1.53t). Alignment angle in respect to the direction of shear was calculated and visualised using Roseplot from R.

Electron microscopy

Liver lobes were perfused and incubated in the fixative solution (2% formaldehyde, 2% glutaraldehyde, 1 mM CaCl2, 1 mM MgCl2, 100 mM Na-phosphate, pH of 7.2) for 2 h at 4 °C. Vibratome sections of 180 μm nominal thickness were post-fixed in 1% osmium tetroxide (buffered in 40 mM Na-Phosphate), dehydrated through graded steps of ethanol and embedded in epoxide (SERVA). Ultrathin sections of 50 nm nominal thickness were post-stained with uranyl and lead for imaging with a Zeiss EM 910 (Zeiss).

Single cell RNA-sequencing

Single cell RNA-sequencing of FACS sorted hepatocyte and LSEC from Double-iECKO and control mice was performed using the CEL-seq2 protocol52. In brief, cells were sorted into 96-well plates containing 2.3 μl of lysis buffer composed of dNTPs (Invitrogen), 1% TritonX-100, ultra-pure water and barcoded primers, centrifuged at maximum speed and snap frozen on dry ice. For hepatocytes, 14 plates containing 1344 cells from four different Double-iECKO mice and other 14 plates consisting of 1344 cells from three different littermate control mice were prepared. For LSEC, 1344 cells from three different mutant mice and 1344 cells from three different control mice were prepared. Cells were lysed and released RNA was annealed to barcoded primers for 3 min at 65 °C. First strand synthesis was initiated using a reaction buffer consisting of 0.4 μl of First strand buffer (Invitrogen), 0.2 μl of DTT (Invitrogen), 0.1 μl of RNase inhibitor (Invitrogen), and 0.1 μl of Superscript IV (Invitrogen) and samples were incubated for 1 h at 50 °C followed by 10 min of heat inactivation at 80 °C. Subsequently, second-strand synthesis was performed according to the original CEL-Seq2 protocol. Barcoded cDNA from each well was pooled per plate and purified using AMPure XP Beads (Beckman Coulter). For linear amplification, the cDNA was subjected to in vitro transcription using MEGAscript T7 (Invitrogen) for 13 h at 37 °C. Excessive primers were removed using EXO-SAP (Thermo Fisher Scientific) and amplified RNA was fragmented using fragmentation reagents (Invitrogen). Amplified RNA fragments were purified using RNAClean XP Beads (Beckman Coulter) and sequencing libraries were prepared according to the original protocol. To account for the low RNA content of LSEC, LSEC libraries were amplified cDNA using 18 cycles for the final PCR, whereas hepatocyte libraries were amplified using 14 cycles, respectively. Both female and male mice were included in the experiment.

Sequencing library was purified from excess primers using AMPure XP Beads and quality controlled using 2100 Bioanalyzer (Agilent) and Qubit fluorometer (Invitrogen). Multiplexed libraries containing 14 pooled libraries were sequenced on a NextSeq550 (Illumina) using Paired-End 75 bp High Output kit to a depth of approximately 100,000 reads per cell. For sequencing, Read1 was set to a length of 13 bases, whereas Read2 was set to 50 bases. 5% PhiX (Illumina) was added to increase complexity and for error rate control.

Single-cell RNA-seq analysis

Raw FASTQ files were demultiplexed by merging R1 reads containing the cell barcode and UMI with their corresponding R2 reads. Subread (v1.6.5)88 was employed to align the reads to the mm10 mouse reference genome. Next, the data was trimmed using Samtools (v1.0.9)89 and count matrices were generated by counting the transcripts for each gene based on the number of UMIs.

Quality control and data analysis were performed using the Seurat V3 package90. In brief, cells were filtered based on the number of detected genes (hepatocytes > 1000 genes and <7000 genes; LSEC > 500 genes and <4000 Genes) and log-normalised. Next, to account for batch effects, data integration was performed. For this, variable genes were detected in the dataset using the FindVariableFeatures function with default settings. A set of anchor genes for correction between batches was computed using the FindIntegrationAnchors function, and biological replicates were integrated using the IntegrateData function. The integrated dataset was scaled using ScaleData followed by principle component analysis using RunPCA. Dimension reduction was performed with the RunUMAP function using 10 PCA dimensions. Clustering of the data was performed using the FindNeighbors with 10 PCs and the FindClusters function with the resolution parameter set to 0.15. Differentially expressed genes (DEG) between the clusters were computed by Wilcoxon rank sum test using the FindMarkers function. Cell clusters that contain immune cell markers and stellate cell markers based on their DEGs were manually removed. Cluster markers were used for cluster annotation based on the marker genes that are enriched for the distinct hepatic zones (Supplementary Fig. 6f)17. Vlnplot and Dotplot functions were employed for the visualisation of the candidate genes.

For the density plots, gene scores were computed using the AddModuleScore function. Density plots of the gene scores of individual cells were generated using ggplot2 (v 3.3.4) and Violin plots were plotted using the Vlnpot function in Seurat. The gene lists for computing gene scores were assembled from previously published literature16,17,63,70,91.

The central hepatocyte gene score was computed based on the expression of Cyp2c29, Mup11, Cldn2, Ang, Gstm1, Inmt, Lect2, Gulo, Cyp2e1, Cyp1a2, Aldh1a1 and Gstm4.

The midzonal hepatocyte gene score was computed based on the expression of Scd1, Mt1, Gm26870, Gm17535, Gm10801, Neu1, Igfbp2 and Cyp8b1.

The portal hepatocyte gene score was computed based on the expression of Hal, Ftcd, Hsd17b6, Pck1, Hsd17b13, Aqp8, Arl4d, Elavl4, Etnppl, Alb, Cyp2f2, Arg1 and Sfxn1.

The liver detoxification gene score was computed based on the expression of Cyp2c37, Gstm2, Cyp2c29, Cyp2e1, Gstm6, Gstm3, Cyp4a12b, Ccyp2d9, Cyp2c37, Gstm2, Fmo1, Gstt1, Gsta3, Cyp2d40, Gclm, Ugt3a2, Cyp8b1, Cyp2a12 and Cyp2f2.

The central LSEC gene score was computed based on the expression of Wnt2, Wnt9b, Rspo3, Kit, Thbd, Slco2a1, Gas6, Igals1, Cd81 and Ptgs1.

The portal LSEC gene score was computed based on the expression of Dll4, Angpt2, Efnb2, Itbp4, Ly6a, Msr1, Efnb1, Ntn4, Galnt15 and Glul.

The shear gene score was computed based on the expression of Itpr1, Sytl2, Cebpd, Eln, Fgf18, Mturn, Fnip2, Fbln2, Neo1, Dusp1, Itga1, Ramp2, Ntn1, Aqp1, Lyst, Srl, Slc22a4, Mmp24, Rapgef4, Plcb4, Kcnn3, Jph4, Hmga2, Fen1, Fjx1, Uhrf1, Hmga1, Map4k4, Dkk1, Slc1a1, Sae1, Apln, Arhgap18, Arhgap24, Junb, Arhgef9, Cxcl3, Hmcn2, Arhgap21, Arhgap23, Arhgef7, Klf2, Klf4, Jun and Fos.

The gap junctional molecule gene score was computed based on the expression of Gja1, Gja4, Gja5, Gjb1 and Gjb2.

DEGs between genotypes were identified using DESeq292 with a two-tailed Wald test and Benjamini–Hochberg. For this, pseudo-bulks for both hepatocytes and LSEC were formed, reflecting the biological replicates and DEG analysis was performed. DEGs were considered significant for adjusted p < 0.05 and log2 fold change > 1.0. DEGs were also calculated in between the spatial bin of the datasets and the respective spatial bins between the genotypes. For this, cells were classified based on the zonation score (see detailed description above) and pseudobulks reflecting the central, midzonal and portal hepatic zones were formed.

Gene set enrichment analyses (GSEA)93,94 was performed on DEGs ranked by fold change using the clusterProfiler v4.0.1 package in _R_95. Gene sets were considered enriched for normalised enrichment scores > 1.5, p.val < 0.05 and q.val < 0.1 for hepatocyte and q.val < 0.40 for LSEC.

For validation, a publicly available dataset was utilised that was deposited on the gene expression omnibus under the accession code GSE19274255. Cell count matrix of CD45− murine was processed accordingly (Supplementary Fig. 7k).

Assignment of zonation position

To each sequenced cell c_c, we assigned a “zonation position”\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$,{\eta }{c}$$\end{document}ηc, which is a quantity indicating the spatial position of the cell along the central-to-portal zonation axis, as inferred from the cell’s expression of zonation markers. The expression of the portal and central markers was summed for each cell, yielding sums \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${S}{c}^{C}$$\end{document}ScC and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${S}{c}^{p}$$\end{document}Scp, and transformed into a zonation position η_η, such that \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\eta }{c}=0$$\end{document}ηc=0 being the portal extreme, and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\eta }_{c}=1$$\end{document}ηc=1 the central end. Values in between indicate mixed expression and hence an intermediate position. Details on how the expression sums were calculated are given below.

For the hepatocytes, we used a modified marker list which is composed of previously established landmark gene sets17 with six additional genes marked by an asterisk:

Markers for central hepatocytes: \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${L}_{C}$$\end{document}LC = (Alad, Aldh1a1, C6, Nat8f2, Cpox, Csad, Cyb5a, Cyp1a2, Cyp2c37, Cyp2c50, Cyp2e1*, Cyp3a11, Glu*, Gstm1, Hpd, Lect2, Mgst1, Oat, Pon1, Prodh, Rgn, Slc16a10).

Markers for portal hepatocytes: \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${L}_{P}$$\end{document}LP = (Alb*, Apof, Arg1, Apom, Asgr2, Asl*, Ass1*, Atp5a1, Cps1, C1s1, C8b, Cpt2, Cyp2f2*, Tkfc, Eef1b2, Elovl2, Fads1, Fbp1, Gc, Gnmt, Hsd17b13, Ifitm3, Igf1, Igfals, Ndufb10, Pck1, Pigr, S100a1, Serpina1c, Serpina1e, Serpind1, Serpinf1, Trf, Uqcrh, Vtn}.

We set \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${k}{{gc}}$$\end{document}kgc for the UMI count of gene g_g in cell c_c, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${s}{c}={\sum }{g}{k}{{gc}}$$\end{document}sc=∑gkgc for the total UMI count of cell c_c, and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${x}{{gc}}$$\end{document}xgc = \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{k}}{{gc}}/{s}{c}$$\end{document}kgc/sc for the fraction of UMIs of cell c_c for gene g_g. For all marker genes, these fractions were normalised by dividing them by the maximum expression of the gene, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${y}{{gc}}={x}{{gc}}/{x}{g}^{\max }$$\end{document}ygc=xgc/xgmax. The maximum \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${x}{g}^{\max }$$\end{document}xgmax is here taken over all hepatocytes of the same genotype.

For each cell c_c the sum of the maximum-normalised expression of the portal markers, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${S}{c}^{P}={\sum }{g\epsilon {L}{P}}{y}{{gc}}$$\end{document}ScP=∑gϵLPygc, and of the central markers, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${S}{c}^{C}={\sum }{g\epsilon {L}{C}}{y}{{gc}}$$\end{document}ScC=∑gϵLCygc was calculated and a fraction, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\eta }{c}={S}{c}^{C}/({S}{c}^{P}+{S}_{c}^{C})$$\end{document}ηc=ScC/(ScP+ScC) was computed. In contrast to Seurat’s AddModuleScore function, we used a non-lognormalised score.

For LSEC, the same method was utilised using the following gene marker lists:

Markers for portal vein LSEC: Ang, Cd81, Emcn, Gas6, Gja1, Gja4, Id4, Ier3, Lfng, Plvap, Ptgs1, Rspo3, Tcim, Wnt2, Wnt9b.

Markers for central vein LSEC: Angpt2, Cxcl9, Dll4, Efnb1, Efnb2, Esm1, Id2, Insr, Jak1, Lama4, Ltbp4, Ly6a, Msr1, Ntn4, Osmr.

Ranked position score

Given Evi and Rspo3 double knockout influences the hepatic zonation pattern, specifically, it reduced the contrast of expression differences between periportal and pericentral cells, some of the zonal marker expressions were lost in Double-iECKO mice. Yet, other marker genes were still informative, thereby allowing us to assign the zonation position. However, the lost markers diluted to signal by reducing the distribution of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\eta }_{c}$$\end{document}ηc in the Double-iECKO sample.

Furthermore, we assumed that _η_η is monotonously related to the actual position of a cell along the central-to-portal axis. However, there was no reason to assume this relation to be linear, or that it is unaffected by treatment or batch effects. Hence, it seemed prudent to use the _η_η values only to order the cells but not to make use of their absolute values.

Therefore, each cell’s η_η value was replaced with the cell rank \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\eta }{Q}$$\end{document}ηQ, scaled to the unit interval, using the cume_dist function in R. Given a list of _n_n cells with their _η_η values, we ordered the cells by increasing _η_η, and then assigned to each cell a rank score, which is simply the rank (zero-based) in the list, divided by _n-1_n−1, such that the first (most portal-like) cell is scored as 0 and the last (most central-like) cell is scored as 1, while the score increases by _1/n_1/n for each cell in between. We referred this score as etaq.

Importantly, the computation was done separately for each experimental batch. For instance, each cell’s etaq values were ranked within a list comprising cells of the same genotype, from the same mouse, and those processed on the same plate. In this manner, we could avoid technical batch effects influencing the result and even out the signal dilution due to the gene knockout described above.

Smoothing gene expression along zonation

Spline regression

To present the expression of a gene along the zonation axis, we used spline regression utilising generalised linear models of the negative-binomial family with logarithmic link with a 4-dimensional B-spline basis. For spline regression, we used a basis of 4 cubic B splines, as produced in R by the bs (basis spline) function of the R core package splines. Yet, the actual regression was performed by glm.nb, the GLM fit function for negative-binomial GLMs from the MASS package96.

The regression yielded coefficients \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\beta }{{rm}}$$\end{document}βrm, one for each of the 4 spline bases (r=1,..,4)(r=1,..,4) and each mouse m_m. A smoothing spline for the gene of interest was obtained for each mouse m_m as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\mu }{m}({\eta }{Q})=,\exp ({\sum }{r}{\beta }{{rm}}{B}{r}({\eta }{Q}))$$\end{document}μm(ηQ)=exp(∑rβrmBr(ηQ)). Here, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\eta }{Q}$$\end{document}ηQ reflected the x-axis for the plot, running from 0 to 1 and indicating the rank-normalised zonation position, and the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${B}{r}(x)$$\end{document}Br(x) were the four basis functions of the standard 4th degree B spline basis. The fitted value \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mu ({\eta }{Q})$$\end{document}μ(ηQ) was the expected expression fraction of the gene in a cell at the position rank score \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\eta }_{Q}$$\end{document}ηQ.

For visualisation, the y-axis depicting \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mu ({\eta }_{Q})$$\end{document}μ(ηQ) was shown with a square-root scaling to improve homo-scedasticity as the square root is the variance-stabilising transformation for the Poisson distribution.

A spatial map of vascular Wnt can be accessed at https://papagei.bioquant.uni-heidelberg.de/kiapp/.

Quantification and statistical analysis

Statistical analysis and graphics were performed with GraphPad Prism v8 and Rstudio v1.2.5042. Unpaired two-tailed Student’s t-tests were performed when comparing two groups. In cases where more than two groups were being compared, one-way ANOVA followed by Bonferroni test was used. Statistical details are described in each figure legend. For pseudobulk DEG and its downsteam GSEA analysis, the samples were processed with a two-tailed Wald test with Benjamini–Hochberg correction. For the proteomics dataset, the samples were tested using a two-tailed t-test, followed by a Benjamini–Hochberg correction. Data are either presented as mean ± standard deviation (s.d.) or mean ± standard error (s.e.m.).

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.