Single-cell RNA-sequencing separates epidermal and mesophyll cell types in the petunia petal.

(A) Representative petunia wild-type (WT), star, wico and phdef flowers from an upper (top picture) and side view (bottom picture). The limb and tube are indicated. Scale bar = 1 cm. Below the flower pictures, a schematic petal primordium is depicted with L1 (future epidermis) and L2 (future mesophyll) cell layers, and PhDEF expression in the WT, star, wico and phdef primordia is indicated in orange. (B) Overview of the experimental protocol to produce petal protoplasts ready for isolation with the 10X Genomics Chromium device. Upper picture: a WT flower cut into small fragments in the cell wall-digesting solution (scale bar = 1 cm). Bottom picture: isolated protoplasts from a WT petal after a 5h-long digestion process, viewed by bright-light microscopy (scale bar = 50 µm). (C) Uniform Manifold Approximation and Projection (UMAP) plot of 11,632 WT petal cells sequenced for their transcriptome, after integration from two biological replicates. Clusters of cells are color-coded and are ordered from the biggest to the smallest, based on the number of cells. (D) Dotplot of key marker genes for petal identity, stamen and carpel identity, epidermal identity and pigmentation. The size of the dot represents the percentage of cells that express a given gene, and the color scale indicates the average expression level of the gene across all cells in a given cluster. st. = stamen, ca. = carpel. (E) Barplot of the number of genes enriched in the epidermis (magenta) or in the mesophyll (green), as defined by varying the log2FoldChange (log2FC) threshold and the layer specificity threshold (% of cells expressing the gene in the epidermis - % of cells expressing the gene in the mesophyll). The Kolmogorov-Smirnov two-sided test (KS test) was applied to compare the gene number distributions between epidermis- and mesophyll-enriched genes at a given layer-specificity threshold, n.s. non significant, * p < 0.05, ** p < 0.005, *** p < 0.001.

PhDEF regulates a different set of genes in the petal epidermis and mesophyll.

(A) Uniform Manifold Approximation and Projection (UMAP) plots of 11,632 WT, 3,875 star and 3,737 wico petal cells sequenced for their transcriptome, after integration and clustering. Here, epidermal clusters (purple) and mesophyll clusters (green, without the vasculature) were merged. (B) Percentage of cells from the different epidermal and mesophyll identities identified from WT, star and wico petal scRNA-Seq. (C) Per-cell expression levels of PhDEF and PhGLO1 in WT, star and wico displayed on UMAP plots. (D) Dotplot of PhDEF and PhGLO1 expression level (color- coded) and percentage of cells expressing the gene (coded in the size of the dot) in the epidermis and mesophyll cells of WT, star and wico petals. (E) Pearson’s correlation coefficient (r) plot between the pseudo-bulk transcriptomes from the mesophyll, epidermis and vasculature cells extracted from the WT, star and wico scRNA-Seq data. A Fisher r-to-z transformation test indicates that all correlation coefficients are significantly different (p < 0.05). (F) Venn diagram of the number of epidermis (purple) and mesophyll (green) DEGs, revealing the number of epidermis- and mesophyll-specific DEGs, as well as common DEGs (i.e. differentially expressed in both layers), based on WT, star and wico scRNA-Seq data. The percentage of activated (light grey) and repressed (dark grey) genes is also displayed. (G) Differential expression (log2FC = log2(FoldChange)) of the common DEGs in the epidermis or in the mesophyll, showing that almost all DEGs are either activated or repressed in both layers. Three genes of interest (PhDEF, PhGLO1 and PhGLO2) are displayed as color points. (H) Ten most enriched Gene Ontology (GO) terms for biological processes in epidermis-specific (left), mesophyll-specific (middle) and common (right) DEGs, after GO term redundancy reduction and sorting by p-value.

PhDEF binds to more target genes in the epidermis than in the mesophyll.

(A) Metaplot of PhDEF binding sites in WT, wico and star petals, with the distance from the Transcription Start Site (TSS) and the Transcription Termination Site (TTS). (B) Heatmap of the read coverage of all peaks detected in WT samples, and of corresponding regions of the genome in wico and star samples. The peaks are sorted according to the highest to lowest read coverage in WT #2, so that the same regions of the genome are in the same line in all librairies. The read coverage is color-coded, and position 0 represents the start of the peak. (C) Read coverage for PhDEF binding to the genomic regions of PhDEF and PhGLO1 in WT, star and wico petals. Peaks that pass our pipeline for reproducibility (see Methods) are displayed as thick red lines. Predicted MADS-binding sites (MADS-bs) are indicated by blue lines, while no TCP-bs were found in these sequences. Although the peak in the promoter region of PhDEF does not pass the IDR threshold, its read coverage profile strongly suggests that PhDEF binds to this position. (D) Intersection of the number of gene-associated peaks from WT, star and wico ChIP-Seq data, and definition of the possible PhDEF binding profiles. The intersection occasionally resulted in the artifical duplication of peaks, and once corrected, this marginally changed the total number of peaks identified for each genotype in the intersection. (E) Distribution of PhDEF binding profiles across its target genes. (F-G) PhDEF binding profile over the genomic locus of AN1 (F) and EOBII (G). Predicted MADS-binding sites (MADS-bs) are indicated by blue lines, and predicted TCP binding sites (TCP-bs) are indicated by orange lines.

Intersection between PhDEF binding and regulatory profiles, and motif enrichment under PhDEF binding sites.

(A) Distribution of PhDEF binding profiles across genes differentially expressed (DEGs) in star only, wico only, or both in star and wico (common DEGs), as compared to WT, from bulk RNA-Seq of petals at stage 8. Stars indicate a significant deviation from expected PhDEF binding profiles, either for the whole distribution (red star, Chi2 goodness-of-fit test) or an enrichment of individual categories (black star, test of equal proportions, one-sided). Binding categories with less than 5 genes were not tested for significant enrichment. (B) Distribution of PhDEF binding profiles across epidermis-specific, mesophyll-specific or common (in both layers) DEGs, obtained from WT, star and wico scRNA-Seq by comparing layer-specific transcriptomes. (C) Selected motifs enriched under PhDEF ChIP-Seq peaks for WT, wico and star petals, bs = binding site. The full list of motifs detected is in Figure S9.