An LKB1–SIK Axis Suppresses Lung Tumor Growth and Controls Differentiation

Recurrently altered tumor suppressors represent key nodes in the programs that constrain tumorigenesis (1). Uncovering the molecular context in which these genes function can aid in understanding how these alterations derail normal cellular physiology and could help in identifying therapeutic opportunities (2). However, the upstream regulators and downstream effectors of a given tumor suppressor are not always identifiable from human cancer genomics data alone, and molecular profiling along with biochemical studies only reveal interactions that could play a role in tumor suppression. Thus, functional genomics studies are ultimately required to directly interrogate the importance of interactors within networks in which bona fide tumor suppressors operate (3).

LKB1 (also known as serine/threonine kinase 11 or STK11) is frequently inactivated in multiple human cancer types, and germline LKB1 mutations play a causative role in the familial cancer susceptibility disorder Peutz–Jeghers syndrome (4, 5). Notably, in lung adenocarcinoma, Lkb1 is mutationally inactivated in 15% to 30% of cases (6). Lkb1 inactivation dramatically accelerates tumor growth in mouse models of lung cancer and broadens the spectrum of histologic subtypes that emerge (7). In a complex with the pseudokinase STRAD and scaffolding protein MO25, LKB1 phosphorylates and activates members of several subfamilies of AMPK-related kinases (AMPK, BRSK, MARK, NUAK, and SIK) through which it governs multiple cellular processes that could, in principle, contribute to tumor suppression (8, 9).

LKB1 and its substrates have been implicated in a broad range of biological processes, including aging, neuronal differentiation, glucose homeostasis, cell polarity, and adhesion (10–15). However, despite multiple lines of evidence supporting functional parallels between LKB1 and several of its substrates, AMPK is conventionally viewed as the focal point of tumor-suppressor activity downstream of LKB1 (9, 16, 17). If AMPK were indeed the principal tumor suppressor downstream of LKB1, then its ablation in vivo would enhance tumor growth. However, AMPK activity has recently been shown to be required for oncogenic KRAS-driven tumorigenesis in the lung, suggesting that other substrates likely suppress tumor growth downstream of LKB1 (18). Thus, the key tumor-suppressive effectors downstream of LKB1 and the contexts in which they operate remain poorly defined.

Here, through CRISPR/Cas9-mediated combinatorial genome editing within a genetically engineered mouse model of lung adenocarcinoma, we investigate the tumor-suppressive capacity of canonical LKB1 substrate families in vivo. We identify SIK1 and SIK3 as potent tumor suppressors and uncover multiple parallels between SIK family inactivation and LKB1 deficiency that underscore the importance of an LKB1–SIK tumor-suppressive axis.

To select candidate LKB1 substrates for characterization in vivo, we examined their expression in normal murine lung epithelium and purified cancer cells from murine lung tumors at different stages of malignant progression (19). The brain-specific substrate-encoding genes Brsk1 and Brsk2, as well as Mark1, were expressed at very low levels in normal lung as well as lung tumors and were thus excluded [Transcripts per million (TPM) <2; Fig. 1A; Supplementary Fig. S1A]. We also examined the expression of these substrates in human lung adenocarcinoma and normal human lung (20). The relative expression of LKB1 substrates in human lung adenocarcinoma largely agreed with their expression in murine tumors, including the low expression of BRSK1/2, the lower expression of AMPKA2 relative to AMPKA1, the lower expression MARK1 relative to the other MARKs, and the higher expression of SIK1 relative to the other SIKs (Fig. 1B; Supplementary Fig. S1B). Beyond expression, analysis of the frequency and nature of mutations in LKB1 substrates in human lung adenocarcinoma did not highlight any substrates as candidate tumor suppressors (Supplementary Fig. S1C). From these analyses, we proceeded to interrogate the tumor-suppressive capacity of the expressed members of the AMPK, MARK, NUAK, and SIK families.

Genetically engineered mouse models enable the genetics of human cancers to be functionally interrogated within autochthonous tumors that recapitulate many histologic and molecular features of human malignancies (21). The recent integration of CRISPR/Cas9-mediated somatic genome editing into these models has accelerated the rate at which genes can be functionally investigated in vivo (3). Additionally, the capacity of CRISPR/Cas9 to simultaneously target multiple loci should enable the study of gene families wherein there may exist functional redundancy. To inactivate LKB1 substrates in vivo using CRISPR/Cas9, we identified single-guide RNAs (sgRNA) that efficiently target each of the 10 expressed LKB1 substrates (Supplementary Fig. S2A; Supplementary Table S1). To address potential functional redundancy between paralogs, we designed lentiviral vectors encoding up to three sgRNAs that simultaneously target the expressed members of the AMPK, MARK, NUAK, and SIK families. In cell culture, we confirmed that the multi-sgRNA design led to gene-targeting efficiencies comparable to vectors encoding single sgRNAs (Supplementary Fig. S2B–S2D).

To assess the tumor-suppressive function of the AMPK, MARK, NUAK, and SIK families of kinases in Kras^G12D-driven lung tumors in vivo, we initiated tumors in Kras^LSL-G12D/+; R26^LSL-Tomato;H11^LSL-Cas9 (KT;H11^LSL-Cas9) mice with each double or triple Lenti-sgRNA/Cre vector (Fig. 1C). Three months after tumor initiation, negative control KT mice (which lack the H11^LSL-Cas9 allele) transduced with Lenti-sgLkb1/Cre and KT;H11^LSL-Cas9 mice transduced with a dual-sgRNA control vector (Lenti-sgControl/Cre) had very few large surface tumors as detected by fluorescence imaging (Fig. 1D and E). In contrast, positive control KT;H11^LSL-Cas9 mice transduced with Lenti-sgLkb1/Cre developed many large tumors (Fig. 1D and E). Strikingly, KT;H11^LSL-Cas9 mice transduced with Lenti-sgSik1–3/Cre also developed many large tumors, whereas KT;H11^LSL-Cas9 mice transduced with lentiviral vectors targeting the MARK, NUAK, and AMPK families developed fewer large tumors. Quantification by histology confirmed that KT;H11^LSL-Cas9 mice transduced with Lenti-sgSik1–3/Cre had high overall tumor burden and many large tumors, indicating that SIKs are tumor suppressors (Fig. 1F and G; Supplementary Fig. S3A).

To confirm that indels were generated at the targeted loci in vivo, we FACS-isolated Tomato^positive neoplastic cells from large individual tumors and bulk tumor-bearing lungs. PCR amplification and sequencing of the targeted regions confirmed that neoplastic cells from tumors initiated with Lenti-sgRNA/Cre vectors targeting the AMPK, NUAK, and MARK families harbored indels at each targeted locus (Supplementary Fig. S3B–S3D). Using this approach, we also observed high-frequency indels (>70%) in Sik1 and/or Sik3 in large Lenti-sgSik1–3/Cre-initiated tumors, indicating that they may be the predominant tumor suppressors of the SIK family (Fig. 2A). In contrast, indels in Sik2 were infrequent, which agrees with previous reports of protumorigenic functions of SIK2 in other cancer types (22–24). Consistent with indels in Sik1 and Sik3, IHC staining for SIK1 and SIK3 revealed reduced protein expression in tumors initiated by Lenti-sgSik1–3/Cre in KT;H11^LSL-Cas9 mice (Fig. 2B and C; Supplementary Fig. S4A).

To assess the tumor-suppressive activity of each SIK paralog, we initiated lung tumors in KT;H11^LSL-Cas9 mice with lentiviral vectors targeting each SIK family member individually, each pairwise combination of SIKs, as well as all three SIKs (Supplementary Fig. S4B). In agreement with our initial results, Lenti-sgControl/Cre-transduced KT;H11^LSL-Cas9 mice and Lenti-sgSik1–3/Cre-transduced KT mice developed very few large tumors (Fig. 2D and E). KT;H11^LSL-Cas9 mice transduced with vectors encoding sgSik1 or sgSik1 and sgSik3 had an increased number of large surface tumors, greater overall tumor burden, and larger individual tumors (Fig. 2D and E; Supplementary Fig. S4C–S4E). In line with our prior indel analyses, these findings indicated that SIK1 and SIK3 are the predominant tumor suppressors of the SIK family.

Given the fundamental importance of the p53 tumor suppressor in lung cancer and previous data suggesting that SIK1 promotes p53-dependent anoikis in transformed human mammary epithelial cells, we determined whether the tumor-suppressor function of SIKs in lung adenocarcinoma was mediated through p53 (20, 25). To test this, we initiated tumors in KT;Trp53^flox/flox (KPT) and KPT;H11^LSL-Cas9 mice with Lenti-sgControl/Cre and Lenti-sgSik1–3/Cre vectors. Lenti-sgSik1–3/Cre-initiated tumors were significantly larger than tumors in any of the control cohorts, suggesting that SIKs mediate tumor suppression even in the context of p53 deficiency (Fig. 2F–H; Supplementary Fig. S5A). Consistent with the KT;H11^LSL-Cas9 context, SIK1 and SIK3 were frequently inactivated in KPT;H11^LSL-Cas9 mice as determined by indel analysis and IHC (Supplementary Fig. S5B and S5C). To confirm that enhanced tumor growth was driven by on-target effects, we used a dual Lenti-sgRNA/Cre vector encoding a second set of sgRNAs targeting Sik1 and Sik3 (Lenti-sgSik1/3^2nd/Cre) to initiate tumors in KPT;H11^LSL-Cas9 mice. KPT;H11^LSL-Cas9 mice transduced with Lenti-sgSik1/3^2nd/Cre had more large tumors and greater overall tumor burden than KPT mice transduced with the same viral vector (Fig. 2I–K; Supplementary Fig. S5D). Collectively, these findings demonstrate that SIKs, particularly SIK1 and SIK3, have tumor-suppressive function that is at least partially independent of p53 activity in our model.

To enable quantitative and multiplexed analyses of tumor suppressors in vivo, we recently integrated CRISPR/Cas9-based genome editing with tumor barcoding and high-throughput barcode sequencing (Tuba-seq; refs. 26, 27). To measure tumor size, Lenti-sgRNA/Cre vectors are diversified with a vector-specific identifier (sgID) and a random barcode (BC), which stably label each transduced cell and all neoplastic cells within the resulting clonal tumors (Fig. 3A; Supplementary Fig. S6A). High-throughput sequencing of the two-component barcodes amplified from bulk genomic DNA of tumor-bearing lungs enables per-tumor measurements of cancer cell number (Fig. 3B; ref. 26). Analysis of the number of cancer cells in tumors of different genotypes relative to tumors initiated with vectors encoding inert sgRNAs can uncover the magnitude of tumor suppression.

Using Tuba-seq, we directly quantified the extent to which the targeting of Sik1 and/or Sik3 phenocopies LKB1 loss in promoting tumor growth. To uncover potential additive or synergistic interactions, we initiated tumors in KT and KT;H11^LSL-Cas9 mice with a pool of barcoded vectors encoding combinations of sgRNAs (Fig. 3A and B). Consistent with previous work, Lkb1 targeting dramatically increased tumor size (Fig. 3C and D; Supplementary Fig. S6B; refs. 26, 27). Permutation of the position of Lkb1-targeting sgRNAs within the dual- or triple-sgRNA configurations resulted in tumor size distributions that were highly correlated, suggesting that sgRNA activity is not strongly influenced by its position within the Lenti-sgRNA/Cre vectors (Supplementary Fig. S6C and S6D).

Targeting of Sik1 led to a modest increase in tumor size, whereas Sik3 targeting had little to no impact on tumor size (Fig. 3C and D; Supplementary Fig. S6B). Strikingly, concomitant targeting of Sik1 and Sik3 was synergistic, increasing tumor size to a magnitude comparable to Lkb1 targeting (Fig. 3C and D). This result was consistent across a series of combinations of two sets of sgRNAs targeting Sik1 and Sik3, suggesting that SIK1 and SIK3 operate as functionally redundant tumor suppressors (Fig. 3C and D). In contrast, targeting of Sik2, in addition to Sik1 and Sik3, did not significantly increase tumor size, thus reinforcing SIK1 and SIK3 as the predominant tumor-suppressive paralogs (Supplementary Fig. S6E and S6F). Finally, in line with our initial screen, targeting of the other families of LKB1 substrates, including a more divergent member of the AMPK family, Snrk, did not dramatically increase tumor growth (Fig. 3E and F; Supplementary Fig. S6G; ref. 28). Together, our findings establish the tumor-suppressive activity of the SIK family, with SIK1 and SIK3 serving as the primary contributors to tumor suppression.

To understand at the cellular level the driving force underlying the increased tumor burden resulting from SIK inactivation, we quantified proliferation and cell death within SIK-targeted and control tumors. IHC staining for phosphorylated histone H3 (Ser10) and Ki-67 revealed a significant increase in the fraction of proliferating cells within SIK-targeted tumors (Fig. 4A–C). Cleaved caspase-3 and TUNEL staining revealed no dramatic change in cell death (Fig. 4A, D, and E). Thus, increased proliferation in SIK-targeted tumors likely underlies their increase in size relative to control tumors.

Similar to their SIK-proficient counterparts, the majority of Lenti-sgSik1–3/Cre- and Lenti-sgSik1/3/Cre-initiated tumors in KT;H11^LSL-Cas9 and KPT;H11^LSL-Cas9 mice were adenomas and early adenocarcinomas. However, most KPT;H11^LSL-Cas9 mice with tumors in which all Sik paralogs were targeted also harbored one or more large regions of invasive acinar adenocarcinoma (Fig. 4F and G). These areas were characterized by the presence of cancer cell clusters and glandular structures embedded within extensive desmoplastic stroma, resembling tumors in KPT;Lkb1^flox/flox mice and previously reported acinar tumors in mouse models of LKB1-deficient lung cancer (Fig. 4F; refs. 29, 30). Importantly, this histologic pattern was not present in KP or KPT mice, even at very late time points, nor was it present in p53-proficient tumors in KT and KT;Lkb1^flox/floxmice or Sik-targeted tumors in KT;H11^LSL-Cas9 mice (Fig. 4G; Supplementary Fig. S7A). In KP tumors, NKX2-1 and HMGA2 are established markers of well-differentiated and poorly differentiated areas, respectively (19, 31). Regardless of their apparent differentiation state, cancer cells within acinar-structured Sik-targeted and LKB1-deficient tumors exhibited heterogeneous expression of NKX2-1 and HMGA2 (Supplementary Fig. S7B). Thus, in addition to enabling the emergence of acinar adenocarcinoma, this uncoupling of differentiation markers is a shared feature of LKB1-deficient and Sik-targeted tumors.

LKB1 loss has been shown to broaden the histologic spectrum of lung tumors within mouse models, enabling the emergence of tumors with squamous, large-cell, and mucinous differentiation (7, 29, 30, 32). Interestingly, LKB1 mutations are frequent in human mucinous lung adenocarcinoma, and gene-expression signatures of mucinous adenocarcinoma are enriched in multiple cohorts of LKB1-mutant human lung adenocarcinomas (Supplementary Tables S2 and S3; refs. 20, 33–38). Patients with Peutz–Jeghers syndrome have also been documented to develop mucinous lung adenocarcinoma (39, 40). KT;H11^LSL-Cas9 and KPT;H11^LSL-Cas9 mice transduced with Sik-targeting lentiviral vectors frequently developed small regions of mucinous lung adenocarcinoma (Fig. 4H and I; Supplementary Fig. S7C–S7E). The abundant intracellular mucin in these regions stained strongly with Alcian blue and phenocopied mucinous regions in KT;Lkb1^flox/flox and KPT;Lkb1^flox/flox mice as well as in KT;H11^LSL-Cas9 mice with Lenti-sgLkb1/Cre-initiated tumors (Supplementary Fig. S7D and S7E). In contrast, mucinous lesions were not detected in KP, KPT, or KT mice nor within KT;H11^LSL-Cas9 or KPT;H11^LSL-Cas9 mice bearing tumors initiated by control vectors (Fig. 4I; Supplementary Fig. S7C). Histology and IHC for p63 and CK5 indicated that neither LKB1 nor SIK deficiency resulted in squamous differentiation, which is in line with lentiviral Cre vectors preferentially generating adenocarcinomas in Kras^LSL−G12D/+;Lkb1^flox/flox mice (Supplementary Fig. S8A–S8D; ref. 29). Together, these analyses demonstrate that targeting of the SIK family recapitulates the expanded histologic spectrum elicited by LKB1 deficiency in lung cancer, which supports the existence of an LKB1–SIK tumor-suppressive axis.

Given the identification of the SIK family of kinases as tumor-suppressive LKB1 substrates, we investigated the impact of LKB1 and SIK inactivation on global gene-expression programs. To uncover the gene-expression changes induced by either LKB1 deficiency or SIK targeting, we isolated Tomato^positive cancer cells from autochthonous tumors and performed RNA sequencing (RNA-seq). Analysis of the reads aligning to the targeted loci confirmed the presence of indels at Sik1 within all five SIK-targeted tumors (Supplementary Fig. S9A).

Hierarchical clustering of gene-expression profiles revealed clustering of SIK-targeted tumors with LKB1-deficient tumors (Fig. 5A). Principal component analysis also uncovered partial clustering of LKB1-deficient and SIK-targeted tumors (Supplementary Fig. S9B). Initial comparison of gene-expression changes elicited by LKB1 deficiency and SIK targeting demonstrated extensive overlap of differentially expressed genes, with more than 70% of the gene-expression changes in LKB1-deficient tumors being recapitulated in the SIK-targeted tumors (Fig. 5B–D; Supplementary Fig. S9C and S9D). Furthermore, direct comparison of all gene-expression changes relative to KT tumors indicated a highly significant correlation between LKB1-deficient and SIK-targeted tumors (Supplementary Fig. S9E).

From a broader perspective, analysis of dysregulated pathways also demonstrated general agreement between the programs that are altered by LKB1 and SIK deficiency. Consistent with the LKB1–SIK axis constraining proliferation in lung tumors, pathways related to proliferation were upregulated in LKB1-deficient and SIK-targeted tumors (Fig. 5E; refs. 7, 41). However, many shared upregulated pathways were unrelated to proliferation, suggesting functional overlap beyond a general increase in proliferation (Fig. 5E). Analysis of motifs enriched at the promoters of genes that were significantly upregulated in LKB1-deficient and SIK-targeted tumors relative to KT tumors revealed an enrichment of CREB/ATF family binding sites (Supplementary Fig. S9F). Several known CREB target genes were significantly upregulated in both LKB1-deficient and SIK-targeted tumors (Supplementary Fig. S9G; refs. 42–44). These findings align with the well-established roles for LKB1 and SIKs as negative regulators of cAMP-driven transcriptional changes (45).

To evaluate the gene-expression overlap between LKB1 and SIK deficiency using independent data sets, we generated a signature composed of genes that are differentially expressed in SIK-targeted tumors relative to control KT tumors. With this signature, we analyzed published microarray-based gene-expression data sets of lung adenocarcinomas generated from different oncogenic KRAS-driven genetically engineered mouse models (7, 30, 46). Applying gene set enrichment analysis (GSEA), we observed a highly significant enrichment of genes that are higher in SIK-targeted tumors among those genes that are higher in LKB1-deficient tumors relative to LKB1-proficient tumors (Fig. 5F; Supplementary Fig. S10A). The genes that are lower in SIK-targeted tumors were also enriched among those genes that are lower in LKB1-deficient tumors relative to other tumor genotypes (Supplementary Fig. S10A). Overall, LKB1-deficient tumors were more similar to SIK-targeted tumors than tumors of other highly proliferative genotypes such as those with inactivation of Trp53 or p16^Ink4a (Supplementary Fig. S10A). Together, these analyses indicate that the transcriptional changes resulting from SIK targeting are uniquely related to the LKB1-deficient state.

To extend our gene-expression analyses to human lung adenocarcinoma, we performed single-sample GSEA on The Cancer Genome Atlas (TCGA) lung adenocarcinoma cohort to score individual patients on the basis of the genes in our SIK signature that were either higher or lower in the SIK-targeted state (20). Upon stratification of patients on the basis of LKB1 genotype, those genes that were lower in SIK-targeted tumors most appropriately changed in LKB1-mutant tumors relative to all other tumors (Fig. 5G; Supplementary Fig. S10B). Alternatively, classification of LKB1-deficient tumors on the basis of a previously published signature of LKB1 deficiency revealed that this subset of tumors exhibited significantly lower expression of SIK-dependent genes (Supplementary Fig. S10C; ref. 47). Consistent with the derivation of our SIK signature from an oncogenic KRAS-driven lung cancer model, the lower scoring of LKB1-mutant tumors relative to all others was conserved in human tumors with oncogenic KRAS alterations (Supplementary Fig. S10D).

Taking an inverse approach, stratification of the distribution of SIK activity scores across patients into tertiles revealed a significant enrichment of LKB1-mutant tumors among those with the lowest SIK signature scores (Fig. 5H). Notably, the few SIK1-mutant tumors in this cohort were binned within the two lowest SIK signature strata, in agreement with these tumors lacking functional SIK1 activity (Supplementary Fig. S10E). Our signature of SIK-dependent genes was also strongly enriched among genes that were lower in LKB1-mutant lung adenocarcinomas relative to all other tumors in three independent data sets as well as in LKB1-mutant human lung adenocarcinoma cell lines relative to LKB1-proficient cell lines (Supplementary Fig. S10F; refs. 48–52).

We also examined whether SIK-targeted tumors exhibit features that have been previously associated with LKB1 loss in human lung cancer and mouse models thereof. Specifically, we used GSEA to assess the enrichment of gene sets generated from previous gene-expression profiling studies of LKB1 deficiency (Supplementary Tables S4 and S5; refs. 47, 49, 53–55). FOX/CREB-regulated genes were significantly higher in SIK-targeted relative to KT tumors (47). Genes related to mitochondrial function as well as molecular chaperoning and stress response were also modestly enriched in SIK-targeted tumors (49). Additionally, several chemokines and cytokines that have been previously shown to be altered in the LKB1-deficient state, such as Ppbp, Il33, Cxcl5, Lgasl9, Ccl5, and Cxcl12, were also dysregulated in SIK-deficient tumors (Supplementary Table S6; ref. 53). In line with frequent downregulation of PD-L1 expression in LKB1-mutant human lung tumors and cell lines, we also observed lower expression of Pdl1 (Cd274) in SIK-targeted tumors (Supplementary Table S6; refs. 49, 56). More generally, SIK-targeted tumors also exhibit significant enrichment of those genes that are upregulated in KRAS-mutant, LKB1-deficient human tumors and cell lines relative to the KRAS-mutant, p53-deficient state (ref. 54; Supplementary Table S5). Collectively, these findings demonstrate extensive parallels at the molecular level between LKB1-deficient human lung tumors and SIK-targeted tumors, further supporting the operation of SIKs within an LKB1-regulated tumor-suppressive axis.

Beyond examining the enrichment of genes dysregulated by SIK deficiency in LKB1-deficient tumors, we analyzed gene-expression data sets derived from LKB1-deficient human lung adenocarcinoma cell lines in which LKB1 had been reexpressed or LKB1-proficient cell lines in which LKB1 had been knocked down (7, 47, 57). The genes that were lower in the SIK-targeted state were highly enriched among those upregulated upon reexpression of LKB1 in three LKB1-deficient cell lines across two independent studies (Fig. 5I; refs. 7, 47). The genes that were higher in the SIK-targeted state were dramatically enriched among those downregulated upon reintroduction of LKB1 into LKB1-deficient A549 cells (Fig. 5I; ref. 47). Complementary to these findings, the genes that were lower in SIK-targeted tumors were strongly enriched among those upregulated upon LKB1 knockdown in LKB1-proficient lung cancer cells (Fig. 5I; ref. 57). Together, these observations suggest that LKB1 regulates genes that are dysregulated in the SIK-deficient context, thus further reinforcing the functional overlap between LKB1 and SIKs.

Our gene-expression analyses substantiated a functional association between SIKs and LKB1; however, we also found that SIK targeting resulted in more extensive remodeling of the transcriptional landscape than LKB1 loss (Fig. 5C and D; Supplementary Fig. S9C–S9E). To address whether SIKs possess additional tumor-suppressive function beyond that controlled by LKB1, we determined whether SIK inactivation conferred a growth advantage to LKB1-deficient tumors. We initiated tumors in KT;Lkb1^flox/flox;H11^LSL-Cas9 mice with either Lenti-sgControl/Cre or Lenti-sgSik1–3/Cre (Supplementary Fig. S11A). Tumor burden, as assessed by lung weight, surface tumor number, and tumor area, was not significantly different between the cohorts (Supplementary Fig. S11B–S11E). However, the size of tumors in the Sik-targeted setting appeared larger than those initiated with Lenti-sgControl/Cre, suggesting some SIK-mediated tumor suppression in the absence of LKB1 (Supplementary Fig. S11F).

To quantify the effect of SIK targeting in LKB1-deficient tumors, we used Tuba-seq. We initiated tumors in KT;Lkb1^flox/flox;H11^LSL-Cas9 mice using the same pool of dual and triple Lenti-sgRNA/Cre vectors that we initially used in KT;H11^LSL-Cas9 mice (Fig. 3A; Supplementary Fig. S11G). Multiple metrics of tumor size revealed that targeting Sik1 or Sik3 individually did not result in a significant increase. However, three of the four vectors targeting Sik1 and Sik3 as well as the vector targeting all three SIKs increased tumor growth (Supplementary Fig. S11H–S11I). These findings suggest that SIK-mediated tumor suppression is not entirely inactivated by LKB1 loss (Supplementary Fig. S12). Interestingly, consistent with previous reports, we find that Sik1 mRNA levels are generally higher in the LKB1-deficient setting (Supplementary Table S7; refs. 47, 58). This increased SIK1 expression may cooperate with an alternate input upstream of SIKs to sustain some degree of tumor suppression in the absence of LKB1 (Supplementary Fig. S12).

Despite LKB1 being one of the most frequently mutated genes in human lung adenocarcinoma and a well-established familial cancer susceptibility gene, our understanding of the mechanisms by which LKB1 constrains tumorigenesis remains limited. AMPK is often viewed as a principal downstream tumor suppressor given its ability to constrain progrowth processes during periods of energetic deficit (9, 16, 17). However, rather than driving proliferation, the ablation of AMPK renders cells more sensitive to metabolic perturbation and severely impairs subcutaneous tumor growth (17). Furthermore, AMPK has been shown to be required for oncogenic KRAS-driven tumor growth in the lung (18). Given the growing complexity surrounding AMPK function during tumorigenesis, there has emerged a pressing need to determine which LKB1 substrates actually constrain tumor growth in vivo.

To identify LKB1 substrates that suppress tumor growth in vivo, we used CRISPR/Cas9-mediated genome editing in a genetically engineered mouse model of lung adenocarcinoma. By simultaneously targeting multiple paralogs across LKB1 substrate families, we identified the SIKs as potent suppressors of lung tumor growth. Beyond driving tumor growth, we identified extensive parallels between LKB1 and SIK loss at the histologic and gene-expression levels, suggesting that SIKs compose a critical tumor-suppressive arm downstream of LKB1.

Although SIK inactivation drives lung tumor growth in vivo, members of this family are rarely inactivated in human cancers (Supplementary Fig. S1C and data not shown). The functional redundancy between SIK1 and SIK3 may render them less favorable targets of direct genetic inactivation relative to LKB1. Consistent with this notion, obligate components of the LKB1 activation complex, namely STRAD and MO25, also exist as multiple paralogs and are only rarely mutated in human cancer (Supplementary Fig. S1C; ref. 9). Beyond genetic inactivation, SIK activity can also be suppressed at the post-translational level (59). Active cAMP signaling can drive PKA-mediated inactivation of SIKs, and this mechanism of SIK inactivation has recently been shown to promote the growth of pancreatic organoids in vitro and after subcutaneous transplantation (60).

Given the position of LKB1 upstream of multiple potentially tumor-suppressive substrates, the direct genetic inactivation of LKB1 likely leads to a greater loss of tumor-suppressor activity than the inactivation of any of its substrates individually. Notably, several LKB1 substrates have been implicated in lung carcinogenesis in other experimental contexts (Fig. 1D–F; refs. 61, 62). It is likely that different substrates play important roles during distinct phases of carcinogenesis; thus, LKB1 inactivation may drive neoplastic cells past multiple barriers during tumor growth and progression in vivo. Finally, other LKB1-regulated substrates outside of the canonical targets investigated here could also possess tumor-suppressive properties.

Beyond defining the compendium of LKB1 substrates, it remains to be fully understood how each of its substrates is regulated and under what circumstances additional pathways contribute to their regulation. Interestingly, our data are consistent with there being residual SIK-mediated tumor suppression in the absence of LKB1, suggesting more complex regulation of SIK function (Supplementary Fig. S12). In line with our findings, alternate activation signals upstream of LKB1 substrates have been described. For instance, AMPK retains some activity in the absence of LKB1 via activation by CAMKK2 under certain conditions. Furthermore, the existence of three kinases upstream of the AMPK ortholog SNF1 in yeast reinforces the potential for multiple regulators of SIKs and other AMPK-related kinases in different settings, including cancer (63, 64).

Another outstanding question concerns what properties drive the difference in tumor-suppressive capacity among the three SIK paralogs. Our data suggest that SIK1 and SIK3, but not SIK2, contribute to LKB1-mediated tumor growth suppression. The three SIK paralogs are structurally quite similar, harboring a highly conserved kinase domain and exhibiting greater divergence at C-terminal regions (59). In contrast to SIK2 and SIK3, which are expressed nearly ubiquitously at a relatively constant level, SIK1 expression is variable and is influenced by processes such as salt intake, hormonal signaling, cellular depolarization, and circadian rhythms (59). Furthermore, due to the limited characterization of the SIK family, the extent to which there exist paralog-specific substrate preferences remains poorly understood.

In line with the findings reported by Hollstein and colleagues in this issue, our results highlight the tumor-suppressive function of a subset of the SIK family of kinases in lung adenocarcinoma (65). Whether these LKB1 substrates suppress the growth of diverse tumors in patients with Peutz–Jeghers syndrome and/or other sporadic cancer types should be investigated. Furthermore, more detailed investigation of the molecular context in which the LKB1–SIK axis of tumor suppression operates may reveal new strategies to pharmacologically counteract or exploit molecular changes that arise from the inactivation of this pathway. Beyond its role in cancer, LKB1 also influences a broad range of normal processes (10–15, 66, 67). Future combinatorial inactivation of SIKs in these biological contexts may uncover a more extensive degree of functional overlap between LKB1 and SIKs.

Lenti-sgRNA/Cre vectors encoding individual sgRNAs were generated as previously described (68). Three to seven sgRNAs were screened for each of the 10 targeted genes. Cutting efficiency was assessed in LSL-YFP;Cas9 reporter cells as described previously (ref. 26; Supplementary Table S1). The primers used for targeted sequencing of targeted loci are listed in Supplementary Table S1. The final sgRNA sequences cloned into viral vectors are listed in Supplementary Fig. S2A. The Neo2, Lkb1, and nontargeting (NT) sgRNA sequences have been previously described (26, 68).

Multi-sgRNA vectors were generated by inserting additional U6-sgRNA cassettes at the Sbf1/PmeI and PacI/NotI cloning sites flanking the original U6-sgRNA cassettes of existing single-sgRNA vectors. The final vectors and the strategies used to construct them are listed in Supplementary Table S1.

For the multi-sgRNA vectors used for Tuba-seq, individual sgRNAs were instead initially cloned downstream of the U6 promoters in pMJ114 (bovine U6), pMJ179 (mouse U6), pMJ117 (human U6; gifts from Jonathan Weissman: Addgene plasmid # 85995, 85996, and 85997; http://n2t.net/addgene:85995, http://n2t.net/addgene:85996, http://n2t.net/addgene:85997; RRIDs:Addgene_85995, Addgene_85996, Addgene_85997) by site-directed mutagenesis (69). Prior to sgRNA cloning, the loxP site-containing region of the mouse U6 promoter of pMJ179 was reverted to its wild-type sequence. For triple-sgRNA vectors, the BamHI site within the homology region between the bU6 and mU6 modules was destroyed to prevent interference with downstream cloning of the sgID-BC cassette. The U6-sgRNA cassettes to be placed in tandem were then PCR-amplified to append homology arms for Gibson assembly, and then they were inserted in tandem at the site of the original U6-sgRNA cassette within our previously described Lenti-sgRNA/Cre vector (backbone previously linearized by PCR; ref. 68). Primer sets used for cloning individual sgRNAs, linearizing the destination vector, and amplifying the U6-sgRNA cassettes for use in Gibson assembly are summarized in Supplementary Table S1.

The assembled multi-sgRNA vectors were then diversified via addition of sgID-BC cassettes as described previously (26). In brief, unique sgID-BC inserts flanked by the BamHI and BspEI sites were produced via PCR with Lenti-sgRNA/Cre as a template using unique forward primers encoding the sgID-BC region and a universal reverse primer. The sgID-BC amplicons were then digested with BamHI and BspEI and ligated into the Lenti-sgRNA/Cre backbones that were previously linearized using BamHI and XmaI. The resultant colonies were then pooled for each vector prior to plasmid DNA extraction.

Lenti-sgRNA/Cre virus was produced as previously described by cotransfection with packaging vectors (delta8.2 and VSV-G) into 293T cells using polyethyleneimine (26). Viral supernatant was collected at 48 and 72 hours after transfection, concentrated by ultra-centrifugation at 2.50 × 10⁴ RPM for 90 minutes, resuspended in PBS, and titered as described previously (26). Vectors used in pooled experiments were packaged individually to prevent the formation of chimeric species (70).

The positional control vector with sgLkb1^2nd driven by mU6 in the triple-sgRNA configuration (sgNeo1-sgLkb1^2nd-sgNT) was excluded from analyses as it was later found to be recombined, impairing sgLkb1^2nd expression. However, the positional control vector with sgLkb1^1st driven by mU6 in the triple-sgRNA configuration (sgNeo1–sgLkb1^1st–sgNT) functioned as expected.

Kras^LSL−G12D, Trp53^flox, Lkb1^flox, H11^LSL−Cas9, and Rosa26^{LSL−tdTomato} mice have been described (68, 71–74). Lung tumors were initiated by intratracheal delivery of Lenti-sgRNA/Cre vectors (75). For our initial screen for tumor-suppressive substrates, cohorts of KT and KT;H11^LSL-Cas9 mice (2–4 mice/group) were transduced with 7.7 × 10⁴ infectious units (IFU)/mouse of each single- and multi-sgRNA Lenti-sgRNA/Cre vector. Tumors were also initiated in KPT and KPT;H11^LSL-Cas9 mice using 5 × 10³ IFU of each multi-sgRNA Lenti-sgRNA/Cre vector (Fig. 2F–K; Supplementary Fig. S5A–S5D). For generating lung tumors for RNA-seq, tumors were also initiated in KT and KT;Lkb1^flox/flox with 2.00 × 10⁵ IFU Lenti-Cre and 1.00 × 10³ IFU Lenti-Cre, respectively. For Tuba-seq, KT, KT;H11^LSL-Cas9, and KT;Lkb1^flox/flox;H11^LSL-Cas9 mice were transduced with 1.2 × 10⁵, 6 × 10⁴, and 2.4 × 10⁴ IFU pooled Lenti-sgRNA/Cre, respectively. Finally, for assessing SIK-mediated tumor suppression in the absence of LKB1 in a nonmultiplexed format, KT;Lkb1^flox/flox;H11^LSL-Cas9 mice were transduced with 6 × 10⁴ IFU Lenti-sgRNA/Cre. The Stanford Institute of Medicine Animal Care and Use Committee approved all animal studies and procedures.

Lung lobes were fixed in 4% formalin and paraffin embedded. Hematoxylin and eosin, Alcian blue, and Masson's trichrome staining was performed using standard methods. Total tumor burden (tumor area/total area × 100%) and individual tumor sizes were calculated using ImageJ. IHC was performed on 4-μm sections with the Vectastain ABC Kit (Vector Laboratories) using the Cadenza system. The following primary antibodies were used: SIK1 (Santa Cruz Biotechnology: sc-83754), SIK3 (Abcam, ab110987), NKX2-1 (Abcam, ab76013), HMGA2 (Biocheck, 59170AP), RFP/Tomato (Rockland, 600401379), cleaved caspase-3 (Cell Signaling Technologies, 9661S), histone H3 phosphorylated at serine 10 (Cell Signaling Technologies, 9701S), cytokeratin 5 (Abcam, ab52635), p63 (Cell Signaling Technologies, 13109S), and Ki-67 (BD Biosciences, 550609). Sections were developed with DAB and counterstained with hematoxylin. The fractions of H3P- and Ki-67–positive nuclei were quantified using ImageJ, and cleaved caspase-3– and TUNEL-positive cancer cells were quantified by direct counting.

Tumors were dissociated using collagenase IV, dispase, and trypsin at 37°C for 30 minutes as previously described (19). Cells were stained with DAPI and antibodies against CD45 (30-F11), CD31 (390), F4/80 (BM8), and Ter119 (all from BioLegend) to exclude hematopoietic and endothelial cells. FACSAria sorters (BD Biosciences) were used for cell sorting.

Indel analysis of ex vivo FACS-isolated neoplastic cells was performed as follows: gDNA was isolated from at least 1.00 × 10⁴ FACS-sorted Tomato^positive cancer cells using either the DNeasy Blood and Tissue Kit or the AllPrep DNA/RNA/Protein Kit (Qiagen). PCR primers were designed to amplify sgRNA-targeted loci, resulting in 350 to 700 bp amplicons specific to each locus. Amplicons were treated with ExoSAP-IT (Thermo Fisher) and sequenced by Sanger sequencing. Cutting efficiency was determined by TIDE analysis (http://tide.deskgen.com/).

Genomic DNA was isolated from bulk tumor–bearing lung from each mouse following the addition of three spike-in controls (5.00 × 10⁵ cells per control) to enable absolute quantification of cell number using Tuba-seq as described previously (26). Libraries were prepared by single-step amplification of the sgID-BC region from a total of 32 μg of genomic DNA per mouse across eight 100-μL reactions using the Q5 High-Fidelity 2 × Master Mix (New England Biolabs, M0494X). To enable computational removal of chimeric reads that result from index hopping during ultra-deep sequencing, the sgID-BCs were amplified using defined dual-indexing primer pairs with unique i5 and i7 indices. The PCR products were subjected to double-sided purification using Agencourt AMPure XP beads (Beckman Coulter, A63881). Purified libraries were then assessed using the Agilent High-Sensitivity DNA Kit (Agilent Technologies, 5067-4626) on the Agilent 2100 Bioanalyzer (Agilent Technologies, G2939BA). Individual libraries were pooled in a weighted format on the basis of lung weight, and the final pool was cleaned up using Agencourt AMPure XP beads. Libraries were sequenced on an Illumina HiSeq 2500 to generate 150-bp paired-end reads.

Only those reads containing complete sgID-BC cassettes were retained, and read pairs exhibiting mismatches in this region were discarded to minimize the impact of sequencing error. Reads were then aggregated on the basis of their random barcode (clonal identifier), representing individual tumors. Barcodes containing indels were discarded to avoid potential alignment errors and miscalculation of distances between barcodes. For any given barcode pileup, smaller barcode pileups within a distance of two nucleotides were removed to minimize the impact of PCR and sequencing errors. Measures of absolute cell number for each tumor were then calculated by multiplying the read counts for each barcode pileup by the size of the spike-in controls (5.00 × 10⁵ cells) and subsequently dividing by the average number of reads for the three barcodes corresponding to the three spike-in controls.

Statistical analyses were performed exclusively on tumors of greater than 500 cells. Multiple metrics of tumor-size distribution were examined, including various percentiles as well as the maximum-likelihood estimate of the mean assuming a log-normal distribution of tumor size (26). P values and standard deviations were calculated by bootstrapping individual mouse tumors for 1,000 permutations. False discovery rates were calculated using the Benjamini–Hochberg procedure.

Total RNA was prepared from FACS-sorted neoplastic cells ranging from 2.50 × 10⁴ to 6.00 × 10⁴ cells. RNA quality was assessed using the RNA 6000 Pico Assay Kit on the Agilent 2100 Bioanalyzer (Agilent). RNA used for RNA-seq had a mean RIN of 8.6. Total RNA (30 ng/sample) was used for cDNA synthesis using the Ovation RNA-seq System (NuGEN Technologies, Inc.). Two micrograms of NuGEN-amplified double-stranded cDNA was sheared using a Covaris sonicator to an average length of 400 bp and subjected to library preparation using the Illumina TruSeq DNA sample preparation kit according to the manufacturer's protocol and sequenced (2 × 100 bp) on an Illumina HiSeq 4000.

Paired-end RNA-seq reads were aligned to the mm10 mouse genome using STAR under standard input parameters (76). To examine for indels at targeted loci in KT;H11^LSL-Cas9 Lenti-sgSik/Cre tumors, aligned reads were visualized in the Integrative Genomics Viewer (77). The number of reads uniquely aligned to exons of individual genes was determined using “intersection-nonempty” mode with HTSeq against the UCSC KnownGene (mm10) transcriptome (78). The differentially expressed genes between different tumor genotypes were called by DESeq2 using the HTSeq-derived counts as input (79). The significance of overlap in terms of differentially expressed genes was computed using hypergeometric tests. The DESeq2-calculated fold changes were used to generate ranked gene lists for input into GSEA (80). The differentially expressed genes with absolute log₂ fold changes greater than 1 and a false discovery rate less than 0.05 were compiled into gene lists that were input into HOMER and analyzed with the findmotifs.pl to perform motif enrichment (window of −400 bp to +100 bp from the TSS and a motif length up to 12 bp; ref. 81). The differentially expressed genes with absolute log₂ fold changes greater than 2 and a false discovery rate less than 0.05 in the comparison of KT;H11^LSL-Cas9 Lenti-sgSik/Cre and KT tumors were compiled into a signature reflecting the SIK-deficient state that was utilized in the analysis of external gene-expression data sets. For direct comparisons of differentially expressed genes with human tumors, a differential expression cutoff of absolute log₂ fold changes greater than 1 and a false discovery rate less than 0.05 were used.

Gene-expression data derived from tumors in humans and genetically engineered mice as well as human cancer cell lines under accession numbers GSE6135, GSE51266, GSE69747, GSE32863, GSE26939, GSE21581, GSE69552, GSE72094, GSE61913, and GSE63882 were acquired from the National Center for Biotechnology Information Gene Expression Omnibus using the GEOquery package (7, 30, 46–52, 56). Differential expression was computed using limma, and the resulting log₂ fold changes were used to generate ranked gene lists for input into GSEA (80, 82). For direct comparisons of differentially expressed genes, a differential expression cutoff of absolute log₂ fold changes greater than 1 and a false discovery rate less than 0.05 were used.

Scaled estimates of gene-expression data, as well as mutation and copy-number data, from the TCGA lung adenocarcinoma cohort were obtained from the Genomic Data Commons data portal (gdc-portal.nci.nih.gov; ref. 20). Only patients for whom gene-expression data were available and LKB1 genotype could be determined were included in subsequent analyses. Patient tumors harboring homozygous deletions of or non–“Putative Passenger” alterations in LKB1 were regarded as LKB1 mutant. Patient tumors harboring amplifications of or non–“Putative Passenger” alterations in KRAS were regarded as KRAS mutant. Single-sample GSEA of Gene Set Variation analysis was used to calculate scores for individual patient samples based on the collective expression levels of genes within the “lower” or “higher” subsets of the SIK deficiency signature (83). The Kolmogorov–Smirnov test was used to determine significant changes signature score distributions. The χ² test was used to determine significant changes in the representation of LKB1 mutants across strata defined by SIK signature scores.

Next-generation sequencing data for the Tuba-seq and RNA-seq experiments have been deposited at the Gene Expression Omnibus under accession number GSE133896.

I.P. Winters is a co-founder of D2G Oncology, Inc., and has ownership interest (including stock, patents, etc.) in the same. M.M. Winslow has received honoraria from the speakers bureaus of Genentech, Merck, and Amgen, has ownership interest (including stock, patents, etc.) in D2G Oncology, and is a consultant/advisory board member for D2G Oncology. No potential conflicts of interest were disclosed by the other authors.

Conception and design: C.W. Murray, J.J. Brady, M.M. Winslow

Development of methodology: C.W. Murray, C. Li, I.P. Winters

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): C.W. Murray, J.J. Brady, M.K. Tsai, I.P. Winters, R. Tang, L. Andrejka, R.K. Ma

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): C.W. Murray, J.J. Brady, M.K. Tsai, C. Li, R.K. Ma, C.A. Kunder

Writing, review, and/or revision of the manuscript: C.W. Murray, M.K. Tsai, C. Li, I.P. Winters, L. Andrejka, R.K. Ma, C.A. Kunder, M.M. Winslow

Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): M.K. Tsai, R.K. Ma, P. Chu

Study supervision: M.M. Winslow

We thank the Stanford Shared FACS facility and Animal Histology Services for technical support; A. Orantes for administrative support; P. Hollstein and R. Shaw for discussions; and D. Petrov, D. Feldser, L. Attardi, S.-H. Chiou, and members of the Winslow laboratory for helpful comments. C.W. Murray and I.P. Winters were supported by the NSF Graduate Research Fellowship Program. C.W. Murray was additionally supported by an Anne T. and Robert M. Bass Stanford Graduate Fellowship. I.P. Winters was additionally supported by NIH F31-CA210627 and NIH T32-HG000044. J.J. Brady was supported by NIH F32-CA189659. C. Li is the Connie and Bob Lurie Fellow of the Damon Runyon Cancer Research Foundation (DRG-2331). R. Tang was supported by a Stanford University School of Medicine Dean's Postdoctoral Fellowship and a TRDRP Postdoctoral fellowship (27FT-0044). This work was supported by NIH R01-CA175336 (to M.M. Winslow), NIH R01-CA207133 (to M.M. Winslow), and NIH R01-CA230919 (to M.M. Winslow) and in part by the Stanford Cancer Institute support grant (NIH P30-CA124435). RNA-seq data were generated on an Illumina HiSeq 4000 that was purchased with funds from the NIH under award number S10OD018220.

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.