RESEARCH ARTICLE

Transcriptome analysis of mouse and human sinoatrial node cellsreveals a conserved genetic programVincent W. W. van Eif1,*, Sonia Stefanovic2,*, Karel van Duijvenboden1, Martijn Bakker1, Vincent Wakker1,Corrie de Gier-de Vries1, Stephane Zaffran2, Arie O. Verkerk1, Bas J. Boukens1 and Vincent M. Christoffels1,‡

ABSTRACTThe rate of contraction of the heart relies on proper developmentand function of the sinoatrial node, which consists of a smallheterogeneous cell population, including Tbx3+ pacemaker cells.Here, we have isolated and characterized the Tbx3+ cells fromTbx3+/Venus knock-in mice. We studied electrophysiologicalparameters during development and found that Venus-labeled cellsare genuine Tbx3+ pacemaker cells. We analyzed the transcriptomesof late fetal FACS-purified Tbx3+ sinoatrial nodal cells andNppb-Katushka+ atrial and ventricular chamber cardiomyocytes,and identified a sinoatrial node-enriched gene program, includingkey nodal transcription factors, BMP signaling and Smoc2, thedisruption of which in mice did not affect heart rhythm. We alsoobtained the transcriptomes of the sinoatrial node region, includingpacemaker and other cell types, and right atrium of human fetuses,and found a gene program including TBX3, SHOX2, ISL1 and HOXfamily members, and BMP and NOTCH signaling componentsconserved between human and mouse. We conclude that aconserved gene program characterizes the sinoatrial node regionand that the Tbx3+/Venus allele provides a reliable tool for visualizingthe sinoatrial node, and studying its development and function.

KEYWORDS: Sinoatrial node,Workingmyocardium, Transcriptionalprofiling, Electrophysiology, BMP-signaling, Human geneexpression

INTRODUCTIONThe sinoatrial node (SAN) is the dominant pacemaker of the heart,located at the junction of the right atrium (RA) and the superiorcaval vein near the crista terminalis. The SAN consists of pacemakercells that are able to depolarize spontaneously, thereby initiating thecardiac cycle. The SAN is innervated by the autonomic nervoussystem, which enables modulation of the beat rate and output.Fibrosis, aging or loss of pacemaker cells may cause SANdysfunction, resulting in bradycardia and sudden death (Adan andCrown, 2003; Dobrzynski et al., 2007). Ultimately, SANdysfunction requires implantation of an electronic pacemaker. Thegeneration of bona fide pacemaker cells by programming stem cellsor reprogramming endogenous non-pacemaker cells has stimulatedmuch interest, as biological pacemakers may overcome the

shortcomings of the currently used electronic devices (Boinket al., 2015; Cingolani et al., 2017; Protze et al., 2017). The ability tocreate biological pacemakers may be improved by the molecular anddevelopmental characterization of the SAN.

Compared with working cardiomyocytes, SAN cells containpoorly developed sarcoplasmic reticula and sarcomeres (Masson-Pevet et al., 1979; Op ‘t Hof, 1988). Furthermore, the electricalactivity of SAN cells differs from working myocardium (WM) as aresult of differences in expression of genes encoding ion-handling(e.g. channels, accessory units) and other proteins (Dobrzynskiet al., 2007; Marionneau et al., 2005; van Kempen et al., 1991). Forexample, the center of the SAN conducts the impulse slowly andexpresses low levels of Gja5 (Cx40), Gja1 (Cx43) and Scn5a(Nav1.5), a set of genes that provide rapid conduction between atrialcardiomyocytes (Dobrzynski et al., 2007; Marionneau et al., 2005).Expression of Gjb6 (Cx30), Gjd3 (Cx30.2) and Gjc1 (Cx45) in theSAN has also been described (Gros et al., 2009; Kreuzberg et al.,2005). Compared with WM cells, SAN cells have a less-negativediastolic membrane potential, a slower action potential upstrokevelocity and do not possess a stable resting membrane potentialbecause of the absence of the inward rectifier K+ current (IK1),which is carried by the Kir2 channel family (Chandler et al., 2009).In addition, SAN cells exhibit the hyperpolarization-activatedcurrent, which is an inward current carried by K+ and Na+ ions(Brown et al., 1979), and is also referred to as ‘funny current’ (If ).This current may contribute to the diastolic depolarization phase ofSAN cell action potentials (DiFrancesco, 2006; Verkerk et al.,2009), and is carried by the hyperpolarization-activated cyclic-nucleotide gated (Hcn) gene family, of which Hcn4 is the mostpredominant isoform found in the SAN (Marionneau et al., 2005;Moosmang et al., 2001; Stieber et al., 2003). Owing to its virtualabsence in the WM and its abundance in the SAN area, Hcn4 hasbecome a potent marker for the mature SAN, although its expressionin the mouse and human embryonic hearts, and duringpathophysiological conditions, is much broader (Liang et al.,2015a; Mommersteeg et al., 2007; Nattel et al., 2008; Sizarov et al.,2011). Mutations in HCN4 in humans are associated withbradycardia (den Hoed et al., 2013; Milano et al., 2014; VerkerkandWilders, 2015). Moreover,Hcn4-deficient mice die in utero dueto severe bradycardia (Stieber et al., 2003). In addition to ionchannels, pacemaker activity is also regulated by a tight coupling ofsarcoplasmic reticulum (SR) Ca2+ cycling molecules with theelectrogenic Na-Ca exchanger, also known as the ‘Ca2+ clock’(Lakatta et al., 2010; Mangoni and Nargeot, 2008).

Previous studies have shown that core transcription factors areinvolved in SAN formation (reviewed by van Eif et al., 2018),including Tbx3 (Hoogaars et al., 2007; Mommersteeg et al., 2007;Wiese et al., 2009), Isl1 (Liang et al., 2015b; Tessadori et al., 2012),Shox2 (Blaschke et al., 2007; Hoffmann et al., 2013; Puskaric et al.,2010) and Tbx18 (Mommersteeg et al., 2010; Wiese et al., 2009).Received 26 October 2018; Accepted 20 March 2019

1Department of Medical Biology, University of Amsterdam, Amsterdam UniversityMedical Centers, Amsterdam 1105 AZ, The Netherlands. 2Aix-Marseille University -INSERM U1251, Marseille Medical Genetics, Marseille 13005, France.*These authors contributed equally to this work

‡Author for correspondence (v.m.christoffels@amc.uva.nl)

V.M.C., 0000-0003-4131-2636

1

© 2019. Published by The Company of Biologists Ltd | Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

Moreover, the key cardiac transcription factor Nkx2-5 is selectivelyrepressed in the developing SAN (Espinoza-Lewis et al., 2009;Mommersteeg et al., 2007; Ye et al., 2015). BMP and WNTsignaling might contribute to SAN development (Bressan et al.,2013; Hashem et al., 2013; Puskaric et al., 2010). Lineage tracingreveals that the definitive SAN develops from Tbx3+ cells in theearly heart tube (Mohan et al., 2018). Tbx3 has been shown toactivate the pacemaker gene program and is able to imposepacemaker function on the atria (Bakker et al., 2012; Hoogaarset al., 2007). Furthermore, Tbx3 suppresses genes associated withatrial WM (Bakker et al., 2012; Hoogaars et al., 2007; Singh et al.,2012). Thus, Tbx3 is a potent regulator and marker of the SAN, andcan be considered as a suitable target gene for functional analysis ofpacemaker tissues.To gain insight into the molecular mechanisms underlying

SAN regulation and function, we set out to obtain the transcriptomeof the SAN pacemaker cells and to measure functional parameters ofthe pacemaker cells throughout development. However, the SANcomprises a heterogeneous cell population, consisting of centralpacemaker cells, peripheral pacemaker cells, fibroblasts, endothelialcells and connective tissue (Dobrzynski et al., 2007; Sanchez-Quintana et al., 2005). Because this heterogeneous cellularcomposition interferes with obtaining uncontaminated pacemakercardiomyocytes for analysis, we generated a mouse model thatexpresses the fluorescent protein Venus under control of the Tbx3locus in the Tbx3+ cells. Patch clamp and optical mappingexperiments on isolated Venus+ and Venus− right atrial cells atembryonic, fetal and adult stages were performed to study theelectrophysiological properties of pacemaker cells duringdevelopment, and to validate the mouse model. To determine thegene transcriptional profile of the SAN, we obtained thetranscriptome of FACS-purified Venus+ SAN cells and Venus−

working cardiomyocytes. In addition, the transcriptome of fetalhuman SANs was assessed and compared with that of mouse, inorder to identify conserved genetic networks driving the pacemakergene program. We identified Smoc2 as a novel conserved highlyspecific SANmarker, and inactivated the gene to assess its function.

RESULTSVenus expression marks Tbx3+ pacemaker cells inTbx3+/Venus miceIn order to visualize the cardiac conduction system, the reportergene encoding the yellow fluorescent protein Venus was insertedinto the start codon of Tbx3 (Fig. 1A and Fig. S1A). Whole-mountin situ hybridization revealed expression of Tbx3 and Venus in thesnout, eye, ear, mammary glands, genital tubercle and limbs ofE10.5 wild-type and Tbx3+/Venus embryos, respectively (Fig. 1A).The pattern of fluorescent Venus signal was highly similar to that ofendogenous Tbx3 (Fig. S1D). Because the Venus-coding sequencewas inserted into the translation start site of Tbx3, expressionof Tbx3 from the modified allele is expected to be disrupted.To validate inactivation, we generated homozygous mutants(Tbx3Venus/Venus) by intercrossing heterozygous Tbx3+/Venus mice.Tbx3Venus/Venus embryos died before developmental stage E12.5and exhibited abnormal forelimbs, a double outlet right ventricleand failed liver development (Fig. S1E), in accordance with thephenotype of Tbx3 homozygous mutants described previously(Bakker et al., 2008; Davenport et al., 2003; Frank et al., 2011;Lüdtke et al., 2009; Mesbah et al., 2008). Tbx3+/Venus mice, incontrast, did not show any abnormalities and were viable and fertile.To investigate whether Venus is expressed selectively in the

Tbx3+ SAN, immunohistochemistry was performed on embryonic,

fetal and adult hearts. Venus was co-labeled with antibodieslabeling the SAN (Tbx3 and Hcn4), atrial working cardiomyocytes(Gja5) and cardiomyocytes in general (cTnI). We found that cTnI+/Hcn4+/Tbx3+ cells co-expressed Venus, which indicates SANpacemaker cardiomyocytes are accurately identified by Venusexpression (Fig. 1B and Fig. S2). No Venus expression was detectedin the Cx40+/Hcn4−/Tbx3− atrial working cardiomyocytes. Withinthe SAN region, Venus was not detected in cTnI−/Hcn4−/Tbx3−

non-cardiomyocytes. However, Venus expression was seen in theTbx3+ cardiac ganglia, which is located at the dorsal side of the atrianear the SAN. These results imply that Venus is expressed in apattern identical to that of Tbx3 and that Tbx3+/Venus mice aresuitable to mark the SAN pacemaker cardiomyocytes.

To gain insight into the structure of the SAN, we performed athree-dimensional reconstruction of the Venus+ area and atria(Fig. 1C). Although the SAN is usually depicted as a comma shapedstructure (Boyett et al., 2000; Wiese et al., 2009), the Venus+ areaappears as a ring structure around the superior caval vein (SCV).The SAN ‘head’ is found at the cranial and right side of the SCV,and is connected with the ‘tail’ running from the sulcus terminalis tothe crista terminalis. In addition, Venus+ myocardium forms aconnection between the left and right atrium and is further visible atthe left side of the SCV in the interatrial septum, connecting theSAN with the atrioventricular conduction axis. The existence ofboth left and right tracts – known as the anterior (left) and posterior(right) intermodal pathway – has been described in mouse, rat andhuman previously (Aoyama et al., 1993; Blom et al., 1999; James,2001; Rentschler et al., 2001). To investigate whether SAN size isaffected by reduced Tbx3 expression in Tbx3+/Venus mice, wedetermined the shape and volume of Hcn4+ SAN regions in ND0wild-type and Tbx3+/Venus mice (both N=3). Comparison revealedthat Tbx3 heterozygosity does not result in altered SAN shape orsize (Fig. 1D).

To study whether the site of initial activation of the atria overlapswith the Venus-expressing (SAN) domain, we performed high-density optical mapping on Tbx3+/Venus tissue preparations (n=14),containing the RA, crista terminalis and intercaval area. Mappingthe endocardial (n=8) and epicardial (n=6) sides of the preparationsrevealed initial pacemaker activity in the center of the Venus+

domain (Fig. 2A). To test sympathetic activity, samples wereexposed to 1 µM isoproterenol. Upon stimulation, beating ratesincreased from 395±27 beats per minute (bpm) to 542±23 bpm(P<0.05), indicating that the sensitivity of the SAN to adrenergicstimulation is maintained in Tbx3+/Venus mice (Fig. S3A).

Patch-clamp experiments were performed on isolated cells toinvestigate the electrophysiological properties. APs of adult Venus+

cells and isolated wild-type SAN cells showed no differences incycle length, RMP, AP duration and DDR, indicating that theheterozygous loss of Tbx3 did not affect the function of the SANcells of Tbx3+/Venus mice (Fig. S2B,C). Venus− cells in the SANregion exhibited APs associated with WM. In contrast, Venus+ cellsshowed spontaneous activity and a diastolic depolarization phase(Fig. 2B). Moreover, Venus+ cells displayed a significantly lessnegative RMP, a lower APA and upstroke velocity, and shorter APdurations at 20% and 50% of repolarization than Venus− cells(Fig. S3B). In addition, we detected an If current in Venus+ cells,whereas a Na+ current (INa; encoded by Scn5a) was not present(Fig. S3C). To study pacemaker cell characteristics duringdevelopment, we compared isolated Venus+ cells from E12.5,E17.5 and adult hearts. At subsequent stages of development, SANpacemaker cells showed increasing automaticity, as displayed by thedecrease in cycle length and increase in DDR (Fig. 2C-D).

2

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

Furthermore, RMP became more negative and APA increased.Finally, APs became shorter and the upstroke velocity increased(Fig. 2C,D).

Transcriptome analysis of the SAN reveals new pacemakermarkers and signaling pathwaysWe obtained transcriptional profiles of both late fetal (E17.5) SANsand working cardiomyocytes using RNA-sequencing. To obtainpurified SAN and working cardiomyocytes, we crossed Tbx3+/Venus

with BAC336-Nppb-Kat mice, in which red fluorescent proteinKatushka is expressed selectively in the cells of the workingcardiomyocytes (Sergeeva et al., 2014). Double transgenic mice(Tbx3+/Venus;BAC336-Nppb-Kat) have Venus+/Katushka− SAN and

Venus−/Katushka+ working cardiomyocytes (Fig. 3A). The SAN,atria and ventricles of double heterozygous E17.5 fetuses wereisolated, dissociated and FACS-purified based on Venus+ orKatushka+ fluorescence (Fig. 3B). Double-positive Venus+/Katushka+ cells were not observed. After performing RNA-seq onthree samples of ±7000 (Venus+) and ±250,000 (Katushka+) cellseach, principal component analysis (PCA) on both tissuepopulations revealed the working cardiomyocyte transcriptomesto be similar, and to be substantially different from the SAN(Fig. S4A). On the other hand, larger variation was observedbetween the SAN samples. The small SAN size and low cellnumbers, along with incomplete cell dissociation before FACSpurification, are possible explanations for sample heterogeneity.

Fig. 1. Venus protein reflects Tbx3 expression inmice. (A) In Tbx3+/Venus mice, a fluorescent reporter(Venus) was inserted in the first exon of the Tbx3-codingregion. Whole-mount in situ hybridization of E10.5 wild-type and Tbx3+/Venus mice, stained for Tbx3 and VenusmRNA, respectively, revealed overlapping expressionpatterns. (B) Immunolabeling of fetal and adult SANregions. Overview images are shown in the first column,merged versions in the last column. Antibodies and colorcodes are as indicated in the images. The boxed areas inthe left-hand images are enlarged in the images on theright. (C) 3D reconstruction of the SAN region and atrialarea. The SAN surrounds the superior caval vein and isformed along a cranial-caudal axis. (D) No difference inSAN volume between wild-type and Tbx3+/Venus micewas found. Data are mean±s.e.m. CGL, cardiacganglion; RA, right atrium; LA, left atrium; SAN, sinoatrialnode; RV, right ventricle; LV, left ventricle; IVS,interventricular septum; SCV, superior caval vein; VV,venus valves; RSCV, right superior caval vein; LSCV, leftsuperior caval vein; PV, pulmonary vein; IT, internodaltract; I, sinoatrial node head; II, sinoatrial node tail orposterior intermodal tract; III, left lateral part of sinoatrialnode; IV, anterior intermodal tract; V, medial intermodaltract; VI, venous valves region. The asterisk indicates thesinus venarum.

3

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

Although distinction between single cells and duplicates is possible,the manual gate setting is arbitrary (Fig. S4B). We identified 16,974genes to be expressed in one or both cell types [normalized reads perkilobase (RPK)>10 reads] (Fig. S4C; Table S1).To verify the purity of the SAN and working cardiomyocyte cell

populations, we compared the expression patterns of SAN andworking cardiomyocyte markers Tbx3, Shox2, Isl1, Gja5, Nppa andNppb between both cells types (Blaschke et al., 2007; Hoffmannet al., 2013; Hoogaars et al., 2007; Liang et al., 2015b; Sun et al.,2006; Vedantham et al., 2015). Read counts for Tbx3 and Shox2were60- and 63-fold enriched, respectively, in theVenus+ cells (Fig. S4C).

HigherGja5, Nppa and Nppb transcripts levels (4-, 15- and 100-foldenriched, respectively) were found in the Katushka+ population(Fig. S4C). Based on the near absence of reads of workingcardiomyocyte markers and very high read counts of SAN markersin each individual SAN transcriptome, we conclude that the quality ofthe SAN samples is adequate (Fig. S4C). We performed clusteranalysis on significantly enriched genes, showing two main groupsseparated by cell type (Fig. 3C). At a significance level ofP<0.05 andlog 2-fold change <−0.5 or >0.5, we identified 3109 SAN-enrichedgenes and 1655 working cardiomyocyte-enriched genes (Fig. 3D,E).Besides the known pacemaker cell markers, we found a subset of

Fig. 2. Electrophysiology of Venus+ cells. (A) Theelectrical impulse originates within the Venus expressiondomain. (B) Typical action potentials recorded from Venus+

and Venus− cells. (C) Typical action potentials of Venus+

cells from E12.5, E17.5 and adult mice. (D) Average actionpotential characteristics of E12.5, E17.5 and adult nodalcells. Data are mean±s.e.m. *P<0.05. RMP, restingmembrane potential; APA, action potential amplitude;APD, action potential duration; DDR, diastolicdepolarization rate.

4

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

Ca2+ channels (Cacna2d2 and Cacna1g), putative Ca2+-bindingproteins (Smoc2 and Vsnl1), ligands and effectors of the BMP-signaling pathway (Bmp2, Bmp4 and Smad9) and neuron-associatedgenes (Th,Cadps andNrp2) to be SAN enriched (Fig. 3D). Functionalannotation on three main clusters revealed that genes expressed inmitochondria (Coq3) and contributing to the contractile apparatus (Ttnand Myh7) and fast conduction (Gja1 and Gja5) were workingcardiomyocyte enriched (Fig. 3C,F). However, cells in the SAN region

expressed higher levels of genes associated with neuronal function(Chrnb4), NOTCH- and BMP-signaling pathways, and calciumhandling (Cacna2d2, Cacna1g) (Fig. 3C,F). Genes functionallyannotated as ‘extracellular region’ and ‘cell adhesion’ were presentin clusters containing working cardiomyocyte and SAN markers,although the latter showed highest significance for these annotationterms (6.5E-19 and 5E-22, respectively). A list of E17.5 SAN- andworking cardiomyocyte-enriched genes is provided in Table S1.

Fig. 3. Transcriptional profiling of E17.5 mouse SAN and WM tissue. (A) Whole-mount fluorescence microscopy of double transgenic E17.5 right atrium.(B) Fluorescence-activated cell-sorting (FACS) report of dissociated nodal and Katushka+ (WM) cells, based onVenus andKatushka fluorescence. (C) Clustering ofgenes, significantly enriched in E17.5 Venus+ (SAN) or WM cells (n=3). Included are genes with over 10 tag counts. Color codes are based on Log2fold changes.SAN-enriched genes appear in blue, genes enriched in WM are depicted in red. (D) Volcano plot, showing gene distributions in the Venus+ population (L2FC<0) orthe WM cardiomyocytes (L2FC>0). (E) 4764 genes (36%) were found to be significantly enriched (P<0.05) in one of both populations. (F) Functional annotationanalysis of SAN- or Chamber-enriched gene clusters depicted in C.

5

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

We used in situ hybridization to validate the expressionpatterns of novel highly SAN region-enriched genes (Alox8,Smoc2, Vsnl1, Gfra2, Chrnb4 and Slc9a3r2) (Fig. 4A).Arachidonate 8-lipoxygenase (Alox8) was found to be expressedpredominantly in the SAN head of E17.5 fetuses, whereas SPARC-related modular calcium binding 2 (Smoc2) and visinin-like protein1 (Vsnl1) marked the entire SAN. Expression of Vsnl1 has beendetected previously in the developing heart at the venous poleregion, although its function remains unknown (Ola et al., 2012).Smoc2 has been shown to be involved in teeth development (Bloch-Zupan et al., 2011) and can be considered as an early intestinal stemcell marker (Muñoz et al., 2012). Interestingly, these novel highlypacemaker-enriched genes have not been associated with SANfunction before. GDNF family receptor α2 (Gfra2) and neuronalacetylcholine receptor subunit β4 (Chrnb4) are expressed in thecardiac ganglia, although Gfra2 also showed expression, to a lesserextent, in the SAN head and tail. Furthermore, we found thatsodium-hydrogen exchange regulatory cofactor 2 (Slc9a3r2) marksthe epithelium in the ventricles and the SAN, including the nodalartery, but to a lesser extend in the atria.The transcriptome analysis indicated BMP signaling to be SAN

enriched. Using in situ hybridization, we found Bmp2 to beexpressed in the Tbx3+/Nppa− SAN region (Fig. 4B,D). Moreover,immunohistochemistry revealed co-staining of phosphorylatedSmad1/5/8 and Tbx3 in the SAN, indicating elevated BMP-signaling activity in the SAN (Fig. 4C). To investigate whetherBMP signaling affects the expression of SAN and WM genes, weincubated dissected fetal RA and LAWM explants and HL-1 cellswith recombinant BMP2 for 48 h and performed qPCR analysis(n=3) (Fig. 3E-G). We found Cacna2d2 expression to besignificantly induced in RA and LA explants after exposure toBMP2 (Fig. 4F). Moreover, incubation of dissected SANs withBmp inhibitor Dmh1 lowered Cacna2d2 levels (Fig. 4F). We alsoanalyzed the expression levels of other SAN andWMmarkers usingRNA sequencing (Fig. 4G). Cacna2d2 and Cacna1g levels weresignificantly elevated, and Bmp2 and Tbx3 was slightly but notsignificantly induced. Conversely, WM markers were slightlyreduced, and Scn5a was significantly decreased after activation ofBMP signaling (Fig. 4G). These observations implicate a role forBMP signaling in activating the pacemaker gene program and thedownregulation of WM genes.Comparing transcriptomes of the SAN derived by independent

preparation methods will improve the accuracy of identifyingtranscriptional profiles of the SAN. Furthermore, analysis of thegene expression profile at different stages provides insight into themechanisms underlying SAN development. Therefore, we comparedour E17.5 samples with the SAN transcriptome published previously(Vedantham et al., 2015) (Table S1). In this study, SANs wereidentified based on GFP expression in Hcn4+ cells and isolated usinglaser capture microscopy (LCM) at different stages (E14.5, P4 andP14). Because the SAN region comprises a mixture of different celltypes (cardiomyocytes, epicardial cells, fibroblasts, endothelial cells,etc.), we expected that these samples would contain non-pacemakercells within the LCM-isolated region. The E14.5 sample in additioncould contain embryonic atrial cells, which still express someHcn4 atthis stage (Liang et al., 2013; Sizarov et al., 2011). Genes that weresignificantly SAN enriched (P<0.05) at one or more stages (E17.5FACS, E14.5 LCM, P4 LCM and P14 LCM), which yielded 6472genes, were clustered. This revealed that genes in clusters associatedwith platelet activation, fibroblast proliferation and neurons areexpressed at higher levels in the E17.5 FACS-purified Venus+ SANpopulation compared with samples isolated by LCM (Fig. 5A).

At all stages, genes involved in cardiac muscle contraction,respiratory chain function and myofibril assembly are moreabundantly expressed in working cardiomyocytes. Effectors of theBMP-signaling pathway (Bmp2 and Bmp4) are SAN enriched inE14.5, E17.5, P4 and P14 mice. The SAN transcriptome is annotatedwith ‘negative regulation of myoblast differentiation’, whereasworking cardiomyocytes express genes associated with positiveregulation of cardiac muscle cell proliferation (Fig. 5A). To study theconsistency between FACS- and LCM-derived samples in moredetail, we assessed the relative expression levels of workingmyocardial, fibroblast, endothelium and neuronal cell markers inthe populations (Fig. S5A-D). Normalized expression levels werecorrected for normalized RPK values of a subset of known andvalidated pacemaker-specific genes (geomean of RPK of Isl1, Tbx3,Shox2, Smoc2, Bmp2 and Bmp4). Although the E14.5 SANsamples (LCM) contained lower levels of transcripts of fibroblast,endothelium and neuronal genes, the working cardiomyocyte markerNppa is present at higher levels compared with other stages,indicating a lower degree of discrimination between SAN andworking cardiomyocyte cells in these samples. Neuronal genes aremore abundantly expressed in the E17.5 FACS Venus+ samplescompared with LCM, which indicates the presence of more neuronalcells. As Tbx3 is expressed in cardiac ganglia (Horsthuis et al., 2009)(Fig. 1B), these cells may have been co-isolated in the Tbx3+/Venus

mouse model.Because the FACS-purified and LCM-isolated samples differ in

cell-type composition, we hypothesize that transcripts consistentlyenriched in both LCM and FACS SAN samples are more likely to bederived from pacemaker cardiomyocytes. We compared expressionprofiles of E14.5 LCM and E17.5 FACS samples, and found 2249genes to be significantly higher expressed in the SAN in bothsamples (>2 fold expression) (Fig. 5B). Functional annotationrevealed that genes involved in Ca2+ handling, NOTCH- andBMP-signaling pathway, and extracellular matrix organization wereSAN enriched (Fig. 5C). In addition, we found a significant numberof genes (>90) associated with neuronal expression. These genesmay be derived from neuronal cells (cardiac ganglia) or could beexpressed by the pacemaker cells themselves (Ebert and Thompson,2001; Gorza et al., 1988; Horsthuis et al., 2009). To investigatewhether the embryonic SAN gene program is maintained at adultstages, we dissected adult SANs and RA (n=2), derived fromTbx3+/Venus;BAC336-Nppb-Kat mice, and found 13,500 genes withRPK-values of at least 10 in the SAN or RAWM tissues (Fig. S6A).Log2-fold changes from dissected adult samples (SAN/WM) werecompared with E14.5 LCM (Fig. S6B) and yielded 379 genes, SANenriched at both stages (L2FC<−1). Tbx3, Isl1, Shox2, Hcn4,Smoc2, Vsnl1 and BMP-signaling components were moreabundantly expressed in the adult SAN, indicating that thesegenes robustly mark the SAN throughout all stages of life.Interestingly, Cntn2 was highly SAN enriched at the adult stagein contrast to the E14.5 stage, at which no expression was detected.We further found Cacna1d to be SAN enriched at both stages, andCacna1b was more abundantly expressed at the adult stage.However, whether Cacna1b levels are derived from neuronalinnervation rather than pacemaker cells remains unknown.

In order to identify transcriptional regulatory componentspossibly involved in SAN development or function, we performedcluster analysis on transcription factors (TFs) of E17.5 FACS,E14.5, P4 and P14 LCM samples that were significantly SAN-enriched at the 1 stage and over (Fig. S7A). For the SAN-enrichedTFs, we found (among others) functional annotations referring tothe BMP- and NOTCH-signaling pathways (Smads), and to cell fate

6

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

commitment (Hoxa2). TFs that were enriched in all datasets(L2FC<0) are listed in Table S1. Besides known TFs (Tbx3, Isl1 andShox2), SAN cells also express Hox-family members (Hoxa1,Hoxa2, Hoxa3, Hoxa4 and Hoxc4), BMP-signaling components(Smad1, Smad6 and Smad9) and Sox genes (Sox8 and Sox10).Pacemaker cell-specific electrophysiological characteristics are

caused by a subset of ion channels in the SAN, the expression ofwhich is different during development, reflecting the alterations inAPs (Fig. 2C,D). To characterize the ion channel distribution in theSAN and WM, we compared expression profiles of ion channel

genes in E17.5 FACS-purified Venus+ and Katushka+ cellpopulations (Fig. S7B). SAN- and chamber-enriched ion channelsare listed in Table S1. We identified 128 SAN-enriched ionchannels (L2FC<0) (Fig. S7B). In contrast to a previous study(Marionneau et al., 2005), we found that Kcna5 (an ultrarapiddelayed rectifier potassium channel), which is activated during theearly repolarization phase, is SAN enriched rather than chamberenriched. As expected, we found Hcn4 and also Hcn1 to be moreabundantly expressed in the SAN. In order to explain the differencesin APs of pacemaker cells (Dobrzynski et al., 2007; Mangoni and

Fig. 4. Transcription analysis reveals new SAN markers and a potential role for BMP signaling in SAN development. (A) In situ hybridization of genesenriched in the Venus+ cells. Chrnb4 marks the cardiac ganglia, whereas Alox8, Smoc2, Vsnl1 and Slc9a3r2 are expressed in the nodal region. Gfra2 wasdetected in both the sinoatrial node and cardiac ganglia. Hcn4 was used as a positive control. Arrowheads indicate regions of expression. (B) In situ hybridizationreveals that Bmp2 is expressed in the Tbx3+/Nppa−SAN region (arrowhead). (C) Phospho-Smad1, -Smad5 and -Smad8 are co-expressed with Tbx3 in the SAN,indicated by arrowheads. (D) Transcriptional profile analysis shows higher levels of Bmp2 in the Tbx3-Venus+ cells compared with the Nppb-Katushka+

population. (E) Strategy to isolate and incubate SAN and atrial myocardium explant tissues with Bmp2 or Dmh1. (F) Real-time RT-PCR analysis revealedincreased expression levels of Cacna2d2 in atrial explants and HL-1 cells after exposure to Bmp2 for 48 h (n=3). Conversely, Dmh1 lowered Cacna2d2 levels inthe SAN tissues, indicating that BMP signaling is involved in the regulation of Cacna2d2 expression (n=3). (G) Transcriptome analysis on atrial explants (n=3),incubated in the presence and absence of Bmp2, shows that BMP-signaling induces transcription of SAN genes and downregulates chamber-enrichedgenes. Cacna1g, Cacna2d2 and Scn5a levels were significantly altered, as indicated by asterisks. Data are mean±s.e.m. in F.

7

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

Nargeot, 2008) during development, we determined the averageexpression levels of ion channels in E12.5, E14.5, P4 and P14LCM samples. Expression of Atp2a2 [which encodes thesarco-endoplasmic reticulum (SR) Ca-ATPase SERCA2A) andPln increased after birth. Expression of Ryr2 and Ryr3, which areresponsible for sarco-endoplasmic reticulum Ca2+ release, wasslightly increased. Expression of Hcn1 and Hcn4 (underlying the Ifcurrent) also increased after birth (Fig. S7B), correlating with thedevelopmental increase of phase 4 depolarization (DDR; Fig. 2).

Expression of Scn5a, which underlies the main cardiac sodiumcurrent, was not altered, whereas Scn10a levels declined. Cacna2d2andCacna1c levels were not notably altered, but we foundCacna1gand Cacnb2 expression to be increased during development.Cacna1g (a T-type Ca2+ channel), in collaboration with Ifcurrents, is responsible for phase 4 depolarization towards the APthreshold. Cacnb2 (a L-type calcium channel) is responsible for theAP upstroke velocity in pacemaker cells, which was increasedduring development (Fig. 2C,D). It is also responsible for an

Fig. 5. Comparison of FACS-purified (E17.5) samples with LCM (E14.5, P4, P14), SAN and WM samples. (A) Heatmap of differentially expressedgenes in E17.5, E14.5, P4 and P14 heart samples with enriched associated GO terms. Color codes are based on Log2fold changes. SAN-enriched genesare shown in blue; WM enriched in red. (B) Gene expression analysis of E17.5 FACS and E14.5 LCM SAN and WM samples. Negative L2FC values depictSAN-enriched genes. Known SAN and WM gene markers are indicated. (C) Functional annotation of SAN-enriched genes.

8

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

important inward current during the AP; thus, the observeddevelopmental AP shortening indicates upregulation ofrepolarizing K+ currents. Indeed, expression of Kcnd3, Kcnh2 andKcnq1 (which underlie transient and rapid and slow delayed rectifierK+ currents) increases during development. We also found anincreased expression of Kcnj3, Kcnj5, Kcnj11, Kcnj12 and Kcnj14,which are responsible for various inward rectifier K+ currents; thisexplains the developmental hyperpolarization of the RMP, despitethe increased expression of Hcn4, Hcn1 and Cacna1g (Fig. 2).

RNA-sequencing of human fetal cells from the SAN regionrevealed a conserved transcription factor programTo identify the transcriptome of human SANs, we isolated 9- and12-week-old (both n=2) human fetal SAN regions and right atrialtissue by microdissection. As in mice, PCA analysis revealed thathuman pacemaker cells show a distinct phenotype compared withmicrodissected non-SAN right atrial tissue (Fig. S8A). We comparedthe log2-fold changes (SAN/RA non-SAN) of 9- and 12-week-oldfetuses and found that 1076 and 698 genes, respectively, were SANregion and RA enriched (L2FC<−0.4) (Fig. 6A). ISL1, SHOX2,TBX3,VSNL1 and SMOC2were enriched in the human SAN regions,as well as genes belonging to signaling pathways (BMP, NOTCH),anterior/posterior pattern specification genes, and genes associatedwith neuronal development and extracellular matrix organization(Fig. 6B). Moreover, expression of genes associated with thecontractile apparatus and fast conduction was lower in the SANarea than in right atrial tissue, in accordance with the expressionprofile in mice. To further compare the SAN region transcriptomes ofhuman and mice, we compared the human and E17.5 mousetranscriptional profiles and identified 855 genes to be SAN region-enriched in both datasets (Fig. 6C). Functional annotation analysisindicates that BMP- and NOTCH-signaling pathway components, aswell as TFs ISL1, SHOX2 and TBX3, and HOX-family members inthe SAN area are conserved between human and mouse. Thepotential role of HOX-family genes is further underlined by theassociation with retinoic acid response (Fig. 6D). We also comparedexpression profiles of E14.5 with human SAN regions, yielding 784enriched genes, and found similar functional annotations (Fig. S8B,C). We found 89 genes to be SAN region enriched at all prenatalstages (mouse E14.5, E17.5 and human 9w+12w; Table S1). Wevalidated the expression levels of SAN region-enriched genes inhuman RA and SAN samples by qPCR (N=4) (Fig. 6F) andconfirmed enrichment of the core program ISL1, SHOX2 and TBX3genes in the SAN as well as members of the BMP family (BMP4and SMAD9), the NOTCH pathway (NOTCH3 and DLL1) andcalcium-binding protein-encoding genes (SMOC2 and VSNL1). Weperformed a gene set enrichment and network analysis (Wang et al.,2017) on mouse E17.5/human fetal SAN region-enriched transcripts(Figs S9, S10 and Table S1). Regulation of development, neuraldifferentiation, extracellular matrix, growth factor binding, calciumand ion channel activity were among the most significantly enrichedbiological processes and molecular functions.

Homozygous Smoc2 mutant hearts do not show alterationsin cardiac electrophysiologyBecause Smoc2 levels are enriched in both murine and human SANs,and its function in the heart unknown, we studied whether inactivationof Smoc2 affects SAN function. UsingCRISPR/Cas9, a 22 bp deletionin exon 3 was made to induce a frameshift mutation. In heterozygousmice, expression levels were reduced by 50% compared with wildtypes (Fig. S11A). In homozygousmice, only lowmRNA levels couldbe detected, probably caused by nonsense-mediated decay

(Frischmeyer and Dietz, 1999). Smoc2−/− mice were found to beviable and fertile, and appeared to be unaffected. We investigatedSmoc2−/−mutants for alterations in cardiac electrophysiology. Both invivo and ex vivo measurements revealed that Smoc2 depletion did notsignificantly affect heart rate (in vivo, both wild type and Smoc2−/−,n=15; ex vivo, wild type, n=7; Smoc2−/−, n=10) (Fig. S11B). Smoc2−/−

mice exhibit a trend towards increased RR intervals (Fig. S11C) orheart rate variability during ex vivo conditions (Fig. S11D), althoughthese differences were not significant. These results indicate that theeffect of Smoc2 inactivation on cardiac function is minimal.

DISCUSSIONWe show that fluorescent reporter Venus targeted to the Tbx3 locusrecapitulates the expression pattern of Tbx3, and can be used toidentify and analyze pacemaker cardiomyocytes. Optical mapping,patch clamp analyses and morphometric analysis indicated that theelectrophysiological properties and size of Venus+ and wild-typeSANs were similar, indicating the presence of the reporter in Tbx3causing heterozygosity is of minor influence. However, theinnervating autonomic nerves also express Tbx3. Therefore, Tbx3heterozygosity may affect neuronal function and, for example, heartrates during stress conditions. Owing to the small size of the SANand its heterogeneous tissue composition, studying the molecularbiology of pacemaker cells has been challenging. EndogenousVenus expression facilitates improved identification, isolation andpurification of pacemaker cells. In recent years, the generation ofbona fide pacemaker cells has gained interest. Generation of abiological pacemaker will benefit from transcriptome analysis ofendogenous SANs. In this study, we used Tbx3+/Venus;BAC336-Nppb-Kat transgenic mice to provide transcriptomes of E17.5 andmature FACS-purified Venus+ SAN and Katushka+ workingcardiomyocytes. Functional annotation clustering was performedto characterize the cell populations. We found expression of genesassociated with anterior/posterior pattern specification, extracellularmatrix organization, neuronal development and BMP- andNOTCH-signaling pathways to be enriched in the SAN. Weidentified novel markers (Smoc2, Vsnl1 and Alox8) not previouslyassociated with SAN development and function, and confirmedtheir expression pattern with in situ hybridization. Inactivation ofSmoc2 by CRISPR/Cas9 genome editing did not result in alteredcardiac function. However, no studies were performed during stressconditions, which could elicit an altered phenotype.

We found that the cluster containing the most highlySAN-enriched genes (Tbx3, Isl1, Shox2, etc.) was functionallyassociated with ‘extracellular region’ and ‘cell adhesion’.Interestingly, functional clustering in a previous study revealedthat genes associated with these processes were downregulated inthe SAN of E12.5 murine embryos when Isl1was inactivated (Lianget al., 2015b). This indicates that genes involved in extracellularmatrix organization and cell adhesion are regulated by thetranscription factor Isl1 in the SAN and may contribute to itsfunction. Moreover, ectopic expression of Tbx2 or Tbx3 in thedeveloping atrial myocardium stimulates BMP and TGFβ signaling,and induces extracellular matrix formation around the Tbx2/3-expressing cardiomyocytes (Singh et al., 2012). This suggests thatT-box factors in the pacemaker cells may function likewise. Theextracellular matrix consists of collagen fibers and forms a regularlypatterned framework in and around the SAN that limits intercellularcontact and limits the conduction velocity (Boyett et al., 2000;Matsuyama et al., 2004; Sanchez-Quintana et al., 2005). Diseases inextracellular matrix organization are associated with cardiacarrhythmias (Adan and Crown, 2003; Thery et al., 1977).

9

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

In contrast to the SAN, WM cells highly expressed genesinvolved in muscle contraction, fast conduction and mitochondrialfunction. Indeed, the SAN contains fewer mitochondria and showsless contractile properties (Boyett et al., 2000; Dobrzynski et al.,2007; Virágh and Challice, 1980). We detected expression levels ofGja1, Gja5 and Scn5a in SAN samples that are not expressed in thenodal center. These transcripts may be derived from tissue withhigher conduction velocity in the SAN periphery, reflecting itsheterogeneity. Therefore, Scn5a+/Gja5+/Hcn4+/Nppa− peripheralcells (Dobrzynski et al., 2005; Lei et al., 2004) may co-express (low)

levels of Tbx3. Furthermore, Scn5a is expressed in the intercavalregion and in larger cells in the SAN rather than in smaller cells(Maier et al., 2003). Single cell sequencing analysis may reveal theextent of heterogeneity of the pacemaker cells, all of which seem toexpress Tbx3 and Shox2. The fraction of Hcn4+ cells that co-expressIsl1 is reduced after birth, indicating that Isl1 expression in the SANbecomes mosaic (Liang et al., 2015b).

We investigated whether the fetal E17.5 murine FACS-purifiedSAN transcriptome is comparable with LCM-derivedtranscriptional profiles of earlier stages of development and

Fig. 6. Transcriptional profile analysis of human fetal SAN regions. (A) Pattern of significantly differentially expressed genes in 9- and 12-week-oldhuman fetuses. Known SAN/RA markers are depicted by arrows. 1076 genes were found to be SAN region enriched at both stages; 698 in the RA (n=2). (B) GOterms and representative genes of SAN- or RA-enriched gene populations. In the SAN region samples, genes associated with neuronal development, NOTCHandBMP signaling are enriched. Genes involved in the contractile apparatus are RA enriched. (C) 855 genes were found to be enriched in both human andmouseE17.5 SANs (L2FC<−0.7). The transcription factors ISL1, SHOX2 and TBX3 are indicated, as well as members of the BMP and NOTCH signaling pathway(BMP2-4, DLL3). (D) GO terms based on genes more abundantly expressed in the SANs of fetal human and E17.5 mice. The NOTCH- and BMP-signalingpathways are represented by SAN cells as well as by genes associated with neuronal development and extracellular matrix organization. (E) 89 genes were foundto be SAN region enriched in bothmurine (E14.5, E17.5) and fetal human (9- and 12-week-old) SAN regions. (F) Real-time RT-PCR validation of RNA-seq resultsobtained with human samples (9 and 12 weeks) (n=4). Data were normalized to HPRT and are expressed as mean fold increase over the working RA samples±s.e.m. *P<0.05.

10

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

adulthood. We found that Isl1, Shox2, Tbx3, Hcn4, Smoc2 andVsnl1, as well as the ion channels Cacna1g, Cacna2d2 and Chrnb4,are SAN enriched from stage E14.5 to adulthood. Interestingly, wefound high levels of Cntn2 in the adult SAN samples, whereas noexpression was detected at E14.5.Cntn2was previously indicated tobe expressed in cardiac conduction system components before birth(Pallante et al., 2010; Qiao et al., 2017). The expression profilingdata implies that Cntn2 expression in the SAN is selective and startsafter birth.Contamination rates, i.e. the presence of non-SAN

cardiomyocytes in the samples, may depend on the techniquesused for cell isolation and purification. Fibroblasts, epicardial cells,endothelial cells, etc. are positioned in close proximity to pacemakercells and might be separated more precisely by FACS purificationthan LCM. However, relatively high levels of fibroblast andendothelial genes were detected in P14 LCM and E17.5 FACSsamples (Fig. S5). This could be explained by the gateway settingprocedure prior to FACS purification, which is arbitrary (Fig. S4B).At E14.5, discrimination between SAN and right atrialcardiomyocytes was less evident due to higher Nppa expressionlevels. Both LCM-isolated and FACS-purified cell populations showabundant expression of neuronal genes, although this was moreprominent in FACS-purified samples. The SAN is innervated byautonomic nervous system cells that play important roles inregulating cardiac function. As Tbx3 was detected in the cardiacganglia (Fig. 1B) (Horsthuis et al., 2009), processing microdissectedand FACS-purified Venus+ SAN cells is more susceptible toneuronal contamination. However, neuronal genes are known to beexpressed in cells of the cardiac conduction system. Neurofilamentwas detected in the SANs of adult mice (Wen and Li, 2015) and inrabbits (Gorza and Vitadello, 1989; Verheijck et al., 1998). Neuronalgene expression was abundant in the atrioventricular canal(precursor of the atrioventricular conduction system) of E10.5murine embryos, long before innervation occurs here (Horsthuiset al., 2009) and in the SAN and bundle branches of rats from E12.5onwards (Nakagawa et al., 1993). Moreover, immunohistochemistryon rabbit heart sections with neurofilament marker iC8 revealedstaining near the atrioventricular junction (Gorza et al., 1988).Furthermore, Pnmt that converts noradrenalin into adrenalin, isexpressed in the SAN of E11.5 rats, days before cardiac innervationoccurs (Ebert and Thompson, 2001). This implies that nervous tissueand developing pacemakermyocardium express an overlapping geneprogram and that the heart is capable of synthesizing and storingcatecholamines itself.Knowledge regarding the molecular biology of pacemaker cells,

which is essential for the generation of biological pacemakers, relieson functional studies in mice and cell lines. Insight into thetranscriptome of human SANs will improve the robustness ofpacemaker cell therapies for application in humans. Owing to thesmall size of the SAN and lack of clear morphological markers,isolation of fetal human SAN cells is arbitrary. Nevertheless, inducedpluripotent stem cells (iPSC) used for generating biologicalpacemakers are immature and therefore transcriptome analysis offetal SAN cells is relevant. We performed transcriptome analysis onmicrodissected human fetal SAN region tissue of 9 and 12 weeks ofgestation, and non-SAN right atrial tissue. We observed that a geneprogram is conserved between human and murine SAN regions,including transcription factors (ISL1, SHOX2, TBX3, Bmp andNotcheffectors and HOXA-family genes) and calcium-binding proteins(SMOC2 and VSNL1). In addition, we found genes associated withneuronal development to be SAN region enriched in both species.Because microdissection of human fetal pacemaker cells is

challenging, a relatively high fraction of atrial and non-pacemakercells was found in the SAN samples, explaining the prominentdifferences in gene transcriptomes between 9- and 12-week-old SANregions (Fig. S7), and the fairly low number of genes consistentlyexpressed at significant higher levels in both SAN area preparations.

We found Notch and Bmp effectors to be SAN region-enriched inhuman andmice. Tbx5, Shox2 and Bmp4were found to play a role ininflow tract development of the mouse embryos (Puskaric et al.,2010) Isl1 deficiency resulted in reduction of Bmp4 expression in theSAN region (Liang et al., 2015b; Vedantham et al., 2015). The role ofBmp4 in SAN development remains elusive. Our data indicatespecific expression of Bmp2/4 ligand and elevated BMP signaling inthe developing SANofmouse and human. CanonicalWNT signalingmight be involved in proliferation of the Tbx18+ mesenchymalprogenitor cell population at the venous pole. The dorsal part of thesinus horns of Ctnnb1-deficient mice was not myocardialized(Norden et al., 2011). However, SAN formation was unaffected(Norden et al., 2011). Canonical WNT signaling was shown to beactive in the SAN and AVC region in E10.5 murine embryos (Gillerset al., 2015). In mouse models, a loss-of-function mutation in WNTsignaling is associated with septum and tricuspid valve defects.Abnormal chamber development and a prolonged PR interval andQRS complex have also been described (Gillers et al., 2015).

Transcription factors belonging to the Hox family were identifiedas SAN region enriched in both murine and human SANs. It hasbeen suggested that initial Hox gene expression is activated byretinoic acid (RA) signaling, as anterior Hox genes exhibit abnormalexpression patterns in embryos with impaired RA synthesis(Raldh2−/− mutant) (Niederreither and Dollé, 2008; Waxman andYelon, 2009). In zebrafish, Hoxb5b acts downstream of RAsignaling and restricts the number of cardiac progenitor cells(Waxman et al., 2008). Moreover, expression patterns of Hoxa1,Hoxa3 and Hoxb1 show overlap with those of Raldh2 in E7.25murine embryos (Bertrand et al., 2011). Mice lacking Raldh2,encoding the enzyme responsible for synthesizing the majority ofembryonic RA, display an impaired expansion of second heart fieldprogenitor cells and morphogenetic defects at the arterial andvenous poles (Ryckebusch et al., 2008). The Hox family genesHoxa1, Hoxa3 and Hoxb1 may be involved in cardiac second heartfield development (Bertrand et al., 2011). First, genetic lineagetracing analysis has shown that Hoxb1+ and Hoxa1+ cells give riseto the inferior wall of the distal OFT and the atrial myocardium(Bertrand et al., 2011). Interestingly, development of the atria andsinus venosus is severely impaired in Raldh2mutants (Niederreitheret al., 2001). Whether this phenotype is caused by depletion of Hoxgene expression has not been fully established.

MATERIALS AND METHODSEthics statementAnimal care was in accordance with national and institutional guidelines.Human fetuses (9 and 12 weeks) were obtained from electively terminatedpregnancies, anonymously donated to research after informed writtenconsent from donors was obtained in concordance with French legislation(PFS14-011) and with prior approval of the protocol (to S.Z.) from the‘Agence de la biomédecine’.

AnimalsGeneration of BAC336-Nppb-Kat mice has been described previously(Sergeeva et al., 2014). A cosmid with Tbx3, isolated from the 129/Olacosmid genomic library obtained from the Resourcenzentrum (RZPD) inBerlin, was kindly provided by Dr Andreas Kispert (Institut fürMolekularbiologie, Medizinische Hochschule Hannover, Germany).Homologous DNA sequences (6.1 kb of upstream sequence and 1.9 kb of

11

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

downstream sequence) were ligated to a Venus-polyA-Frt-flanked PGK-neocassette derived from pKOII (Bardeesy et al., 2002) to generate a Tbx3-targeting construct (Fig. S1A) in which the first three codons of the Tbx3-coding region were replaced by the Venus-pA cassette (Nagai et al., 2002).The linearized targeting construct was electroporated into E141B10embryonic stem (ES) cells to generate targeted cell lines. A diphtheriatoxin A cassette was used to positively select for homologous recombinants.Chimeras were generated by injection of targeted ES cells into C57Bl6 hostblastocysts. Germline transmission of the targeted allele was obtained bymating with FVB females. Subsequently, Tbx3-VenusNEO mice werecrossed with FlpE mice (Rodríguez et al., 2000) to remove the PGK-neocassette. Removal of the Neo cassette was assessed by Southern blottingand PCR (Fig. S1B,C). Progeny were screened by PCR for the presenceof the Tbx3-Venus allele using the following primers: fw1(AGCGGAGCCAAGCCAGCA), rv1 (CCTTGGCCTCCAGGTGCAC)and rv2 (TTGATGCCGTTCTTCTGCTTGT). The Tbx3-Venus allele hasbeen maintained on the FVB genetic background.

CRISPR/Cas9DNA target sites and corresponding oligos were identified using ZiFitTargeter software. Oligos were annealed [95°C for 5 min and 25°C for1 min (ramp down 0.1°C/s)] using 1× T4 ligation buffer (Invitrogen) andligated into BsaI-digested and gel purified sgRNA expression vectorpDR274, using T4 DNA ligase (Invitrogen). DNA was isolated with theJETSTAR DNA isolation kit (ITK, 200050). Expression vectors forsgRNAs and Cas9 (3 µg each) were linearized for 4 h at 37°C with DraI(NEB, R0129S) and PmeI (NEB, R0560S), respectively, followed byphenol-chloroform purification. Transcription of sgRNAs and Cas9 wasperformed using the MEGAshort T7 kit (Life Technologies, AM1354) andmMessage mMachine T7 Ultra kit (Life Technologies, AM1345),respectively, followed by purification with the MEGAclear kit (LifeTechnologies, AM1908), according to the manufacturer’s protocol. For theinjections, we used 5 ng/µl sgRNA and 10 ng/µl Cas9 mRNA.

Cellular electrophysiologyCell isolationSingle cells were isolated from the SAN region and right atria of E12.5,E17.5 or adult hearts by an enzymatic dissociation procedure modifiedfrom Verkerk et al. (2009). In brief, the SAN region and atria were excised,cut into small strips (0.3-0.5×1 mm) and stored in a cold (6°C) modifiedTyrode’s solution containing (in mM): NaCl 140, KCl 5.4, CaCl2 1.8,MgCl2 1.0, glucose 5.5 and HEPES 5.0 (pH 7.4). Next, the strips wereplaced in nominally Ca2+-free Tyrode’s solution (20°C), i.e. modifiedTyrode’s solution with 10 µMCaCl2, which was refreshed twice. Then, thestrips were incubated for 11-13 min in nominally Ca2+-free Tyrode’ssolution (37°C) to which liberase IV (0.25-0.29 U/ml, Roche) and elastase(2.4-0.7 U/ml, Serva) were added. During the incubation period, the stripswere triturated through a pipette (tip diameter 2.0 mm). The dissociationwas stopped by transferring the strips into a modified Kraft-Brühe (KB)solution (37°C), that was refreshed three times (at 20°C). The modified KBsolution contained (in mM): KCl 85, K2HPO4 30, MgSO4 5.0, glucose 20,pyruvic acid 5.0, creatine 5.0, taurine 30, β-hydroxybutyric acid 5.0,succinic acid 5.0, BSA 1% and Na2ATP 2.0 (pH adjusted to 6.9 usingKOH). Subsequently, the strips were stored for 30 min and then triturated(pipette tip diameter: 1.2 mm) in modified KB solution (20°C) for 2 minto obtain single cells. Single cells were stored at room temperature forat least 45 min in modified KB solution before they were put into arecording chamber on the stage of an inverted microscope, and superfusedwith modified Tyrode’s solution (37°C). Individual Venus+, Venus− orwild-type cells were isolated, and individual spindle and elongatedspindle-like cells displaying regular contractions for SAN experimentswere selected.

Data acquisitionMembrane potentials and currents were recorded using the amphotericin-perforated patch-clamp technique using an Axopatch 200B amplifier(Molecular Devices). Signals were low-pass filtered at 10 kHz cut-offfrequency, and digitized at 25 kHz. Data acquisition and analysis were

accomplished using custom software. Pipettes (borosilicate glass;resistance 2–4 MΩ) were filled with a solution containing (in mM): K-gluconate 125, KCl 20, NaCl 5, amphotericin B 0.22 and HEPES 10 (pHadjusted to 7.2 using KOH). Action potentials (APs) were measured in bothspontaneous active and quiescent cells. In quiescent cells, APs were elicitedat 2 Hz by 3 ms current pulses, 1.5× threshold current pulses through thepatch pipette. We analyzed cycle length, maximal diastolic potential(MDP), diastolic depolarization rate (DDR, measured over the 50 ms timeinterval starting at MDP+1 mV), maximal upstroke velocity (Vmax),maximal AP amplitude (APA) and AP duration at 20, 50 and 90%repolarization (APD20, APD50 and APD90). Parameters from 10consecutive APs were averaged.

Optical mappingThe mice were killed by cervical dislocation, after which the heart wasexcised and placed in a solution (30°C) containing (in mM): NaCl 128, KCl4.7, CaCl2 1.45, MgCl2 0.6, NaHCO3 27, NaH2PO4 0.4 and glucose 11 (pHwas maintained at 7.4 by equilibration with a mixture of 95% O2 and 5%CO2). Fat and non-cardiac tissue was removed from the venous pole of theheart. Then the right atrium, crista terminalis and intercaval area weredissected free, whereafter the preparation was pinned down. Subsequently,the isolated preparations were incubated in 10 ml Tyrode’s solutioncontaining 15 μM Di-4 ANEPPS, superfused and placed in the opticalmapping setup. Excitation light was provided by a 5 W power LED (filtered510±20 nm). Fluorescence (filtered>610 nm) was transmitted through atandem lens system on CMOS sensor (100×100 elements, sampling rate5 kHz, MICAM Ultima). Optical action potentials were analyzed usingcustom-made software.

StatisticsData are mean±s.e.m. Two groups were compared using an unpaired t-testand more than two groups were compared using one-way ANOVA followedby a Student-Newman-Keuls method post-hoc test. P<0.05 indicatesstatistical significance.

Electrocardiogram (ECG)Animals were anaesthetized with 5% isoflurane (Pharmachemie) andmaintained under 1.5-2.0%. Electrodes were placed in the right (R) and left(L) armpit, and in the left groin (F). An electrocardiogram (ECG) wasrecorded (PowerLab 26T; AD-Instruments) for a period of 5 min. ECGparameters were determined in Lead II (L-R) based on the last 30 s of therecording. To record ECGs ex vivo, the adult mice were killed by cervicaldislocation. The heart was rapidly excised, cannulated and mounted on aLangendorff perfusion set-up, as described previously (Boukens et al.,2013). Hearts from neonatal mice (ND0-1) were isolated and superfusedwith HEPES-buffered Tyrode’s solution containing (in mM): 140 NaCl, 5.4KCl, 1.8 CaCl2, 1.0 MgCl2, 5.5 glucose and 5.0 HEPES at a temperature of36±0.2°C; pH was adjusted to 7.4 using NaOH. During perfusion, the heartwas submerged and electrodes were placed at the right (R) and left (L) sideof the base of the heart and at the left side of the apex (F) at a 5 mm distance.Lead II was used to determine ECG parameters.

Explant cultures from mouse fetal tissueE17.5 SANs and atrial WM tissues were dissected and placed in M199medium (Sigma-Aldrich, M4530) supplemented with 2% fetal calf serum(Invitrogen, 10270-106), 1% penicillin/streptomycin (Thermo-Fisher,15140-122) and incubated at 37°C. HL-1 cells were cultured at 37°C inClaycomb medium (Sigma-Aldrich, 51800C) containing 1% HL-1 definedFBS (Lonza, LO77227), 1% penicillin/streptomycin (Thermo-Fisher,15140-122), 0.2 mM norepinephrine (Sigma-Aldrich, A0937-1G)(dissolved in 30 mM ascorbic acid solution) and 1% glutamax (Thermo-Fisher, 35050-061). HL-1 cells were seeded on plates coated with 0.02%gelatin and 0.5% fibronectin (Sigma-Aldrich, F-1141). HL-1 cells wererecently authenticated and tested for contamination. Recombinant-HumanBMP2 was added (R&D Systems, 355-BM) to atrial myocardial explantsand HL-1 cells to a final concentration of 200 ng/ml. SAN samples wereincubated with 3 µMDmh1 (Sigma-Aldrich, D8946-5MG). All tissues wereincubated for 48 h and medium was refreshed after 24 h.

12

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

Tissue isolation and fluorescence-activated cell sorting (FACS)for RNA-seqDissected tissue was collected in cold PBS (6°C) and dissociated usingTrypsin/EDTA (0.25%/1 mM) to obtain a single cell suspension. Cells werepassed through a cell strainer to exclude debris and cell clumps.Fluorescence-activated cell sorting (FACS) was performed on a Sony cellsorter (SH800Z, Sony Biotechnology). The human SAN and working atrialmyocardium were recognized by their anatomic landmarks. Under themicroscope, the right atrium was opened to expose the crista terminalis, theintercaval area and the interatrial septum. A thin strip of SAN tissue, limitedby crista terminalis, atrial septum and orifices of the venae cavae was cutfrom the right atrium.

RNA isolation, real-time RT-PCR and RNA-sequencingRNA from isolated SANs, atrial WM cells and cultured explants wasobtained using the Nucleospin RNA XS isolation kit (Bioke, MN740902.50), according to manufacturer’s protocol. The quality of totalRNA was assessed by micro-electrophoresis on acrylamide gel (Agilent2100 Bioanalyser). For real-time RT-PCR Superscript II (Invitrogen) andOligodT primers were used to generate cDNA templates. Expression levelswere assessed with quantitative RT-PCR using the Lightcycler Real-TimePCR system (Roche Diagnostics). Hprt was used as reference gene. ForRNA-seq, the RNA template was processed with the Ovation V2 RNA-seqsystem (Nugen, 7102-32), according to the manufacturer’s protocol.Libraries were generated using the Ovation Ultra-low system V2 1-96(Nugen, 0344NB-A01). Amplified cDNA was bead purified (AmpureXP,Beckman-Coulter). Samples were sequenced on a Hiseq 2500 sequencingsystem (Illumina) with 50 bp single-end reads. RNA-seq data have beendeposited in GEP under accession number available under GSE125932.

ImmunohistochemistryHearts were fixed overnight in 4% formaldehyde, embedded in paraplastand sectioned at 7 μm. Embedding media were removed from sections and,by boiling for 5 min in a high pressure cooking pan in Antigen UnmaskingSolution (Vector H3300), the epitopes in the tissue were unmasked. Tissuesections were incubated with primary antibodies overnight in Tris-NaClbuffer with blocking powder. Primary antibodies used were mouse-anti-troponin 1 (Millipore, MAB1691, 1:400), rabbit-anti-cardiac troponin 1(Hytest LTD, 4T21_2, 1:200), mouse-anti-Cx40 (USBiological, C7856,1:200), goat-anti-GFP (Abcam, 839963, 1:100), rabbit-anti-GFP (SZE2110, 1:100), rabbit-anti-Hcn4 (Chemicon, NG164345, 1:200), goat-anti-Tbx3 (Santa Cruz, B1006, 1:200) and rabbit-anti-phospho-SMAD1/-5/-8 (Cell Signaling, 9511, 1:200). Secondary antibodies used were Alexa680 (1:250; Molecular Probes, A21084), Alexa 568 (1:250; MolecularProbes, A10042) and Alexa 488 (1:250; Molecular Probes, A11017) withthe appropriate epitope to visualize the primary antibody. For detection ofgoat anti-Tbx3 we used biotinylated donkey anti-goat antibody (1:250). Fordetection of biotinylated antibodies, we used the TSA Enhancement kit(Perkin Elmer, NEL702).

Fluorescence microscopyWe used a Leica MZ FL III microscope (filter excitation 470/40 andemission 590LP) and a Leica DFC320 camera to view and photographfluorescent hearts.

Whole mount in situ hybridizationEmbryos were fixed in 4% formaldehyde.Whole mount in situ hybridizationwas performed as described previously (Moorman et al., 2001).

AcknowledgementsWe thank Ingeborg Hooijkaas, Rajiv Mohan and Bouke de Boer for theircontributions.

Competing interestsThe authors declare no competing or financial interests.

Author contributionsConceptualization: V.W.W.v.E., S.S., K.v.D., M.B., A.O.V., B.J.B., V.M.C.;Methodology: V.W.W.v.E., S.S., K.v.D., M.B., V.W., C.d.G.-d.V., S.Z., A.O.V., B.J.B.,

V.C.; Software: K.v.D., M.B., C.d.G.-d.V., A.O.V., B.J.B.; Validation: V.W.W.v.E.,S.S., K.v.D., M.B., A.O.V., B.J.B.; Formal analysis: V.W.W.v.E., S.S., K.v.D., M.B.,A.O.V., B.J.B.; Investigation: V.W.W.v.E., S.S., K.v.D., M.B., V.W., C.d.G.-d.V., S.Z.,A.O.V., B.J.B., V.M.C.; Resources: V.M.C.; Data curation: V.W.W.v.E., S.S., K.v.D.,M.B., V.W., C.d.G.-d.V., S.Z., A.O.V., B.J.B.; Writing - original draft: V.W.W.v.E.,V.M.C.; Writing - review & editing: S.S., M.B., V.W., S.Z., A.O.V., B.J.B., V.M.C.;Visualization: V.W.W.v.E., S.S., K.v.D., M.B., S.Z., A.O.V., B.J.B., V.M.C.;Supervision: V.W., C.d.G.-d.V., S.Z., A.O.V., B.J.B., V.M.C.; Project administration:V.M.C.; Funding acquisition: V.M.C.

FundingS.S. received support from the European Society of Cardiology (ESC training grant2013). B.J.B. received support from the Hartstichting (2016T047). This study wassupported by the Nederlandse Organisatie voor Wetenschappelijk Onderzoek(ZonMW TOP 40-00812-98-12086) and by the Fondation Leducq (14CVD01to V.M.C.).

Data availabilityAll RNA-seq data have been deposited in GEO under accession numberGSE125932.

Supplementary informationSupplementary information available online athttp://dev.biologists.org/lookup/doi/10.1242/dev.173161.supplemental

ReferencesAdan, V. and Crown, L. A. (2003). Diagnosis and treatment of sick sinus syndrome.

Am. Fam. Physician 67, 1725-1732.Aoyama, N., Kikawada, R. and Yamashina, S. (1993). Immunohistochemical

study on the development of the rat heart conduction system using anti-Leu-7antibody. Arch. Histol. Cytol. 56, 303-315. doi:10.1679/aohc.56.303

Bakker, M. L., Boukens, B. J., Mommersteeg, M. T. M., Brons, J. F., Wakker, V.,Moorman, A. F. M. and Christoffels, V. M. (2008). Transcription factor Tbx3 isrequired for the specification of the atrioventricular conduction system. Circ. Res.102, 1340-1349. doi:10.1161/CIRCRESAHA.107.169565

Bakker, M. L., Boink, G. J. J., Boukens, B. J., Verkerk, A. O., van den Boogaard,M., den Haan, A. D., Hoogaars, W. M. H., Buermans, H. P., de Bakker, J. M. T.,Seppen, J. et al. (2012). T-box transcription factor TBX3 reprograms maturecardiac myocytes into pacemaker-like cells. Cardiovasc. Res. 94, 439-449.doi:10.1093/cvr/cvs120

Bardeesy, N., Sinha, M., Hezel, A. F., Signoretti, S., Hathaway, N. A., Sharpless,N. E., Loda, M., Carrasco, D. R. and DePinho, R. A. (2002). Loss of the Lkb1tumour suppressor provokes intestinal polyposis but resistance to transformation.Nature 419, 162-167. doi:10.1038/nature01045

Bertrand, N., Roux, M., Ryckebusch, L., Niederreither, K., Dolle, P., Moon, A.,Capecchi, M. and Zaffran, S. (2011). Hox genes define distinct progenitor sub-domains within the second heart field. Dev. Biol. 353, 266-274. doi:10.1016/j.ydbio.2011.02.029

Blaschke, R. J., Hahurij, N. D., Kuijper, S., Just, S., Wisse, L. J., Deissler, K.,Maxelon, T., Anastassiadis, K., Spitzer, J., Hardt, S. E. et al. (2007). Targetedmutation reveals essential functions of the homeodomain transcription factorShox2 in sinoatrial and pacemaking development. Circulation 115, 1830-1838.doi:10.1161/CIRCULATIONAHA.106.637819

Bloch-Zupan, A., Jamet, X., Etard, C., Laugel, V., Muller, J., Geoffroy, V.,Strauss, J.-P., Pelletier, V., Marion, V., Poch, O. et al. (2011). Homozygositymapping and candidate prioritization identify mutations, missed by whole-exomesequencing, in SMOC2, causing major dental developmental defects.Am. J. Hum. Genet. 89, 773-781. doi:10.1016/j.ajhg.2011.11.002

Blom, N. A., Gittenberger-de Groot, A. C., DeRuiter, M. C., Poelmann, R. E.,Mentink, M. M. T. and Ottenkamp, J. (1999). Development of the cardiacconduction tissue in human embryos using HNK- 1 antigen expression: possiblerelevance for understanding of abnormal atrial automaticity. Circulation 99,800-806. doi:10.1161/01.CIR.99.6.800

Boink, G. J. J., Christoffels, V. M., Robinson, R. B. and Tan, H. L. (2015). Thepast, present, and future of pacemaker therapies. Trends Cardiovasc. Med. 25,661-673. doi:10.1016/j.tcm.2015.02.005

Boukens, B. J., Hoogendijk, M. G., Verkerk, A. O., Linnenbank, A., van Dam, P.,Remme, C.-A., Fiolet, J. W., Opthof, T., Christoffels, V. M. and Coronel, R.(2013). Early repolarization in mice causes overestimation of ventricular activationtime by the QRS duration. Cardiovasc. Res. 97, 182-191. doi:10.1093/cvr/cvs299

Boyett, M. R., Honjo, H. and Kodama, I. (2000). The sinoatrial node, aheterogeneous pacemaker structure. Cardiovasc. Res. 47, 658-687. doi:10.1016/S0008-6363(00)00135-8

Bressan, M., Liu, G. and Mikawa, T. (2013). Early mesodermal cues assign aviancardiac pacemaker fate potential in a tertiary heart field. Science 340, 744-748.doi:10.1126/science.1232877

Brown, H. F., DiFrancesco, D. and Noble, S. J. (1979). How does adrenalineaccelerate the heart? Nature 280, 235-236. doi:10.1038/280235a0

13

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

Chandler, N. J., Greener, I. D., Tellez, J. O., Inada, S., Musa, H., Molenaar, P.,DiFrancesco, D., Baruscotti, M., Longhi, R., Anderson, R. H. et al. (2009).Molecular architecture of the human sinus node: insights into the functionof the cardiac pacemaker. Circulation 119, 1562-1575. doi:10.1161/CIRCULATIONAHA.108.804369

Cingolani, E., Goldhaber, J. I. and Marban, E. (2017). Next-generationpacemakers: from small devices to biological pacemakers. Nat. Rev. Cardiol.15, 139-150. doi:10.1038/nrcardio.2017.165

Davenport, T. G., Jerome-Majewska, L. A. and Papaioannou, V. E. (2003).Mammary gland, limb and yolk sac defects inmice lacking Tbx3, the genemutatedin human ulnar mammary syndrome. Development 130, 2263-2273. doi:10.1242/dev.00431

den Hoed, M., Eijgelsheim, M., Esko, T., Brundel, B. J., Peal, D. S., Evans, D. M.,Nolte, I. M., Segre, A. V., Holm, H., Handsaker, R. E. et al. (2013). Identificationof heart rate-associated loci and their effects on cardiac conduction and rhythmdisorders. Nat. Genet. 45, 621-631. doi:10.1038/ng.2610

DiFrancesco, D. (2006). Serious workings of the funny current. Prog. Biophys. Mol.Biol. 90, 13-25. doi:10.1016/j.pbiomolbio.2005.05.001

Dobrzynski, H., Li, J., Tellez, J., Greener, I. D., Nikolski, V. P., Wright, S. E.,Parson, S. H., Jones, S. A., Lancaster, M. K., Yamamoto, M. et al. (2005).Computer three-dimensional reconstruction of the sinoatrial node. Circulation111, 846-854. doi:10.1161/01.CIR.0000152100.04087.DB

Dobrzynski, H., Boyett, M. R. and Anderson, R. H. (2007). New insights intopacemaker activity: promoting understanding of sick sinus syndrome. Circulation115, 1921-1932. doi:10.1161/CIRCULATIONAHA.106.616011

Ebert, S. N. and Thompson, R. P. (2001). Embryonic epinephrine synthesis in therat heart before innervation: association with pacemaking and conduction tissuedevelopment. Circ. Res. 88, 117-124. doi:10.1161/01.RES.88.1.117

Espinoza-Lewis, R. A., Yu, L., He, F., Liu, H., Tang, R., Shi, J., Sun, X., Martin,J. F., Wang, D., Yang, J. et al. (2009). Shox2 is essential for the differentiation ofcardiac pacemaker cells by repressing Nkx2-5. Dev. Biol. 327, 378-385. doi:10.1016/j.ydbio.2008.12.028

Frank, D. U., Carter, K. L., Thomas, K. R., Burr, R. M., Bakker, M. L., Coetzee,W. A., Tristani-Firouzi, M., Bamshad, M. J., Christoffels, V. M. andMoon, A.M.(2011). Lethal arrhythmias in Tbx3-deficient mice reveal extreme dosagesensitivity of cardiac conduction system function and homeostasis. Proc. Natl.Acad. Sci. USA 109, E154-E163. doi:10.1073/pnas.1115165109

Frischmeyer, P. A. and Dietz, H. C. (1999). Nonsense-mediated mRNA decay inhealth and disease. Hum. Mol. Genet. 8, 1893-1900. doi:10.1093/hmg/8.10.1893

Gillers, B. S., Chiplunkar, A., Aly, H., Valenta, T., Basler, K., Christoffels, V. M.,Efimov, I. R., Boukens, B. J. andRentschler, S. (2015). Canonical wnt signalingregulates atrioventricular junction programming and electrophysiologicalproperties. Circ. Res. 116, 398-406. doi:10.1161/CIRCRESAHA.116.304731

Gorza, L. and Vitadello, M. (1989). Distribution of conduction system fibres in thedeveloping and adult rabbit heart, revealed by an antineurofilament antibody.Circ.Res. 65, 360-369. doi:10.1161/01.RES.65.2.360

Gorza, L., Schiaffino, S. and Vitadello, M. (1988). Heart conduction system: aneural crest derivative? Brain Res. 457, 360-366. doi:10.1016/0006-8993(88)90707-X

Gros, D., Theveniau-Ruissy, M., Bernard, M., Calmels, T., Kober, F., Sohl, G.,Willecke, K., Nargeot, J., Jongsma, H. J. and Mangoni, M. E. (2009). Connexin30 is expressed in the mouse sino-atrial node, and modulates heart rate.Cardiovasc. Res. 85, 45-55. doi:10.1093/cvr/cvp280

Hashem, S. I., Lam, M. L., Mihardja, S. S., White, S. M., Lee, R. J. and Claycomb,W. C. (2013). Shox2 regulates the pacemaker gene program in embryoid bodies.Stem Cells Dev. 22, 2915-2926. doi:10.1089/scd.2013.0123

Hoffmann, S., Berger, I. M., Glaser, A., Bacon, C., Li, L., Gretz, N., Steinbeisser,H., Rottbauer, W., Just, S. and Rappold, G. (2013). Islet1 is a directtranscriptional target of the homeodomain transcription factor Shox2 andrescues the Shox2-mediated bradycardia. Basic Res. Cardiol. 108, 339. doi:10.1007/s00395-013-0339-z

Hoogaars,W.M. H., Engel, A., Brons, J. F., Verkerk, A. O., de Lange, F. J.,Wong,L. Y. E., Bakker, M. L., Clout, D. E., Wakker, V., Barnett, P. et al. (2007). Tbx3controls the sinoatrial node gene program and imposes pacemaker function on theatria. Genes Dev. 21, 1098-1112. doi:10.1101/gad.416007

Horsthuis, T., Buermans, H. P. J., Brons, J. F., Verkerk, A. O., Bakker, M. L.,Wakker, V., Clout, D. E. W., Moorman, A. F. M., ‘t Hoen, P. A. C. andChristoffels, V. M. (2009). Gene expression profiling of the formingatrioventricular node using a novel Tbx3-based node-specific transgenicreporter. Circ. Res. 105, 61-69. doi:10.1161/CIRCRESAHA.108.192443

James, T. N. (2001). The internodal pathways of the human heart. Prog.Cardiovasc. Dis. 43, 495-535. doi:10.1053/pcad.2001.24598

Kreuzberg, M. M., Sohl, G., Kim, J.-S., Verselis, V. K., Willecke, K. andBukauskas, F. F. (2005). Functional properties of mouse connexin30.2expressed in the conduction system of the heart. Circ. Res. 96, 1169-1177.doi:10.1161/01.RES.0000169271.33675.05

Lakatta, E. G., Maltsev, V. A. and Vinogradova, T. M. (2010). A coupled SYSTEMof intracellular Ca2+ clocks and surface membrane voltage clocks controls thetimekeeping mechanism of the heart’s pacemaker. Circ. Res. 106, 659-673.doi:10.1161/CIRCRESAHA.109.206078

Lei, M., Jones, S. A., Liu, J., Lancaster, M. K., Fung, S. S.-M., Dobrzynski, H.,Camelliti, P., Maier, S. K. G., Noble, D. andBoyett, M. R. (2004). Requirement ofneuronal- and cardiac-type sodium channels for murine sinoatrial nodepacemaking. J. Physiol. 559, 835-848. doi:10.1113/jphysiol.2004.068643

Liang, X., Wang, G., Lin, L., Lowe, J., Zhang, Q., Bu, L., Chen, Y., Chen, J.,Sun, Y. andEvans, S.M. (2013). HCN4 dynamicallymarks the first heart field andconduction system precursors. Circ. Res. 113, 399-407. doi:10.1161/CIRCRESAHA.113.301588

Liang, X., Evans, S. M. and Sun, Y. (2015a). Insights into cardiac conductionsystem formation provided by HCN4 expression. Trends Cardiovasc. Med. 25,1-9. doi:10.1016/j.tcm.2014.08.009

Liang, X., Zhang, Q., Cattaneo, P., Zhuang, S., Gong, X., Spann, N. J., Jiang, C.,Cao, X., Zhao, X., Zhang, X. et al. (2015b). Transcription factor ISL1 is essentialfor pacemaker development and function. J. Clin. Invest. 125, 3256-3268. doi:10.1172/JCI68257

Ludtke, T. H.-W., Christoffels, V. M., Petry, M. and Kispert, A. (2009). Tbx3promotes liver bud expansion during mouse development by suppression ofcholangiocyte differentiation. Hepatology 49, 969-978. doi:10.1002/hep.22700

Maier, S. K. G., Westenbroek, R. E., Yamanushi, T. T., Dobrzynski, H., Boyett,M. R., Catterall, W. A. and Scheuer, T. (2003). An unexpected requirement forbrain-type sodium channels for control of heart rate in the mouse sinoatrial node.Proc. Natl. Acad. Sci. USA 100, 3507-3512. doi:10.1073/pnas.2627986100

Mangoni, M. E. and Nargeot, J. (2008). Genesis and regulation of the heartautomaticity. Physiol. Rev. 88, 919-982. doi:10.1152/physrev.00018.2007

Marionneau, C., Couette, B., Liu, J., Li, H., Mangoni, M. E., Nargeot, J., Lei, M.,Escande, D. and Demolombe, S. (2005). Specific pattern of ionic channel geneexpression associated with pacemaker activity in themouse heart. J. Physiol. 562,223-234. doi:10.1113/jphysiol.2004.074047

Masson-Pevet, M., Bleeker, W. K. and Gros, D. (1979). The plasmamembrane ofleading pacemaker cells in the rabbit sinus node. A qualitative and quantitativeultrastructural analysis. Circ. Res. 45, 621-629. doi:10.1161/01.RES.45.5.621

Matsuyama, T. A., Inoue, S., Kobayashi, Y., Sakai, T., Saito, T., Katagiri, T. andOta, H. (2004). Anatomical diversity and age-related histological changes in thehuman right atrial posterolateral wall. Europace 6, 307-315. doi:10.1016/j.eupc.2004.03.011

Mesbah, K., Harrelson, Z., Theveniau-Ruissy, M., Papaioannou, V. E. and Kelly,R. G. (2008). Tbx3 is required for outflow tract development. Circ. Res. 103,743-750. doi:10.1161/CIRCRESAHA.108.172858

Milano, A., Vermeer, A. M. C., Lodder, E. M., Barc, J., Verkerk, A. O., Postma,A. V., van der Bilt, I. A. C., Baars, M. J. H., van Haelst, P. L., Caliskan, K. et al.(2014). HCN4 mutations in multiple families with bradycardia and left ventricularnoncompaction cardiomyopathy. J. Am. Coll. Cardiol. 64, 745-756. doi:10.1016/j.jacc.2014.05.045

Mohan, R. A., Mommersteeg, M. T. M., Domınguez, J. N., Choquet, C., Wakker,V., de Gier-de Vries, C., Boink, G. J. J., Boukens, B. J., Miquerol, L., Verkerk,A. O. et al. (2018). Embryonic Tbx3+ cardiomyocytes form the mature cardiacconduction system by progressive fate restriction.Development 145, 1-12. doi:10.1242/dev.167361

Mommersteeg, M. T. M., Hoogaars, W. M. H., Prall, O. W. J., de Gier-de Vries, C.,Wiese, C., Clout, D. E. W., Papaioannou, V. E., Brown, N. A., Harvey, R. P.,Moorman, A. F. M. et al. (2007). Molecular pathway for the localized formation ofthe sinoatrial node. Circ. Res. 100, 354-362. doi:10.1161/01.RES.0000258019.74591.b3

Mommersteeg, M. T. M., Domınguez, J. N., Wiese, C., Norden, J., de Gier-de,V. C., Burch, J. B. E., Kispert, A., Brown, N. A., Moorman, A. F. M. andChristoffels, V. M. (2010). The sinus venosus progenitors separate and diversifyfrom the first and second heart fields early in development. Cardiovasc. Res. 87,92-101. doi:10.1093/cvr/cvq033

Moorman, A. F. M., Houweling, A. C., de Boer, P. A. J. and Christoffels, V. M.(2001). Sensitive nonradioactive detection of mRNA in tissue sections: novelapplication of the whole-mount in situ hybridization protocol. J. Histochem.Cytochem. 49, 1-8. doi:10.1177/002215540104900101

Moosmang, S., Stieber, J., Zong, X., Biel, M., Hofmann, F. and Ludwig, A. (2001).Cellular expression and functional characterization of four hyperpolarization-activated pacemaker channels in cardiac and neuronal tissues. Eur. J. Biochem.268, 1646-1652. doi:10.1046/j.1432-1327.2001.02036.x

Munoz, J., Stange, D. E., Schepers, A. G., van de Wetering, M., Koo, B.-K.,Itzkovitz, S., Volckmann, R., Kung, K. S., Koster, J., Radulescu, S. et al. (2012).The Lgr5 intestinal stem cell signature: robust expression of proposed quiescent‘+4’ cell markers. EMBO J. 31, 3079-3091. doi:10.1038/emboj.2012.166

Nagai, T., Ibata, K., Park, E. S., Kubota, M., Mikoshiba, K. and Miyawaki, A.(2002). A variant of yellow fluorescent protein with fast and efficient maturation forcell-biological applications. Nat. Biotechnol. 20, 87-90. doi:10.1038/nbt0102-87

Nakagawa, M., Thompson, R. P., Terracio, L. and Borg, T. K. (1993).Developmental anatomy of HNK-1 immunoreactivity in the embryonic rat heart:co-distribution with early conduction tissue. Anat. Embryol. 187, 445-460. doi:10.1007/BF00174420

Nattel, S., Burstein, B. and Dobrev, D. (2008). Atrial remodeling and atrialfibrillation: mechanisms and implications. Circulation. Arrhythm. Electrophysiol. 1,62-73. doi:10.1161/CIRCEP.107.754564

14

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT

Niederreither, K. and Dolle, P. (2008). Retinoic acid in development: towards anintegrated view. Nat. Rev. Genet. 9, 541-553. doi:10.1038/nrg2340

Niederreither, K., Vermot, J., Messaddeq, N., Schuhbaur, B., Chambon, P. andDolle, P. (2001). Embryonic retinoic acid synthesis is essential for heartmorphogenesis in the mouse. Development 128, 1019-1031.

Norden, J., Greulich, F., Rudat, C., Taketo, M. M. and Kispert, A. (2011). Wnt/beta-catenin signaling maintains the mesenchymal precursor pool for murinesinus horn formation. Circ. Res. 109, e42-e50. doi:10.1161/CIRCRESAHA.111.245340

Ola, R., Lefebvre, S., Braunewell, K.-H., Sainio, K. and Sariola, H. (2012). Theexpression of Visinin-like 1 during mouse embryonic development. Gene Expr.Patterns 12, 53-62. doi:10.1016/j.gep.2011.11.004

Op ‘t, and Hof, T. (1988). The mammalian sinoatrial node.Cardiovasc. Drugs. Ther.1, 573-598. doi:10.1007/BF02125744

Pallante, B. A., Giovannone, S., Fang-Yu, L., Zhang, J., Liu, N., Kang, G., Dun,W., Boyden, P. A. and Fishman, G. I. (2010). Contactin-2 expression in thecardiac Purkinje fiber network. Circ. Arrhythm. Electrophysiol. 3, 186-194. doi:10.1161/CIRCEP.109.928820

Protze, S. I., Liu, J., Nussinovitch, U., Ohana, L., Backx, P. H., Gepstein, L. andKeller, G. M. (2017). Sinoatrial node cardiomyocytes derived from humanpluripotent cells function as a biological pacemaker. Nat. Biotechnol. 35, 56-68.doi:10.1038/nbt.3745

Puskaric, S., Schmitteckert, S., Mori, A. D., Glaser, A., Schneider, K. U.,Bruneau, B. G., Blaschke, R. J., Steinbeisser, H. and Rappold, G. (2010).Shox2 mediates Tbx5 activity by regulating Bmp4 in the pacemaker region of thedeveloping heart. Hum. Mol. Genet. 19, 4625-4633. doi:10.1093/hmg/ddq393

Qiao, Y., Lipovsky, C., Hicks, S., Bhatnagar, S., Li, G., Khandekar, A., Guzy, R.,Woo, K. V., Nichols, C. G., Efimov, I. R. et al. (2017). Transient notch activationinduces long-term gene expression changes leading to sick sinus syndrome inmice. Circ. Res. 121, 549-563. doi:10.1161/CIRCRESAHA.116.310396

Rentschler, S., Vaidya, D. M., Tamaddon, H., Degenhardt, K., Sassoon, D.,Morley, G. E., Jalife, J. and Fishman, G. I. (2001). Visualization and functionalcharacterization of the developing murine cardiac conduction system.Development 128, 1785-1792.

Rodrıguez, C. I., Buchholz, F., Galloway, J., Sequerra, R., Kasper, J., Ayala, R.,Stewart, A. F. and Dymecki, S. M. (2000). High-efficiency deleter mice show thatFLPe is an alternative to Cre-loxP. Nat. Genet. 25, 139-140. doi:10.1038/75973

Ryckebusch, L., Wang, Z., Bertrand, N., Lin, S.-C., Chi, X., Schwartz, R., Zaffran,S. and Niederreither, K. (2008). Retinoic acid deficiency alters second heart fieldformation. Proc. Natl. Acad. Sci. USA 105, 2913-2918. doi:10.1073/pnas.0712344105

Sanchez-Quintana, D., Cabrera, J. A., Farre, J., Climent, V., Anderson, R. H. andHo, S. Y. (2005). Sinus node revisited in the era of electroanatomical mapping andcatheter ablation. Heart 91, 189-194. doi:10.1136/hrt.2003.031542

Sergeeva, I. A., Hooijkaas, I. B., van derMade, I., Jong,W.M. C., Creemers, E. E.and Christoffels, V. M. (2014). A transgenic mouse model for the simultaneousmonitoring of ANF and BNP gene activity during heart development and disease.Cardiovasc. Res. 101, 78-86. doi:10.1093/cvr/cvt228

Singh, R., Hoogaars, W. M., Barnett, P., Grieskamp, T., Rana, M. S., Buermans,H., Farin, H. F., Petry, M., Heallen, T., Martin, J. F. et al. (2012). Tbx2 and Tbx3induce atrioventricular myocardial development and endocardial cushionformation. Cell Mol. Life Sci. 69, 1377-1389. doi:10.1007/s00018-011-0884-2

Sizarov, A., Devalla, H. D., Anderson, R. H., Passier, R., Christoffels, V. M. andMoorman, A. F. M. (2011). Molecular analysis of patterning of conduction tissuesin the developing human heart.Circ. Arrhythm. Electrophysiol. 4, 532-542. doi:10.1161/CIRCEP.111.963421

Stieber, J., Herrmann, S., Feil, S., Loster, J., Feil, R., Biel, M., Hofmann, F. andLudwig, A. (2003). The hyperpolarization-activated channel HCN4 is required for

the generation of pacemaker action potentials in the embryonic heart. Proc. Natl.Acad. Sci. USA 100, 15235-15240. doi:10.1073/pnas.2434235100

Sun, Y., Liang, X., Najafi, N., Cass, M., Lin, L., Cai, C.-L., Chen, J. and Evans, S.(2006). Islet 1 is expressed in distinct cardiovascular lineages, includingpacemaker and coronary vascular cells. Dev. Biol. 304, 286-296. doi:10.1016/j.ydbio.2006.12.048

Tessadori, F., van Weerd, J. H., Burkhard, S. B., Verkerk, A. O., de Pater, E.,Boukens, B. J., Vink, A., Christoffels, V. M. and Bakkers, J. (2012).Identification and functional characterization of cardiac pacemaker cells inzebrafish. PLoS ONE 7, e47644. doi:10.1371/journal.pone.0047644

Thery, C., Gosselin, B., Lekieffre, J. and Warembourg, H. (1977). Pathology ofsinoatrial node. Correlations with electrocardiographic findings in 111 patients.Am. Heart J. 93, 735-740. doi:10.1016/S0002-8703(77)80070-7

van Eif, V. W. W., Devalla, H. D., Boink, G. J. J. and Christoffels, V. M. (2018).Transcriptional regulation of the cardiac conduction system. Nat. Rev. Cardiol. 15,1-14. doi:10.1038/s41569-018-0031-y

van Kempen, M. J. A., Fromaget, C., Gros, D., Moorman, A. F. M. and Lamers,W. H. (1991). Spatial distribution of connexin43, the major cardiac gap junctionprotein, in the developing and adult rat heart. Circ. Res. 68, 1638-1651. doi:10.1161/01.RES.68.6.1638

Vedantham, V., Galang, G., Evangelista, M., Deo, R. C. and Srivastava, D.(2015). RNA sequencing of mouse sinoatrial node reveals an upstream regulatoryrole for Islet-1 in cardiac pacemaker cells. Circ. Res. 116, 797-803. doi:10.1161/CIRCRESAHA.116.305913

Verheijck, E. E., Wessels, A., van Ginneken, A. C. G., Bourier, J., Markman,M. W. M., Vermeulen, J. L. M., de Bakker, J. M. T., Lamers, W. H., Opthof, T.and Bouman, L. N. (1998). Distribution of atrial and nodal cells within the rabbitsinoatrial node: models of sinoatrial transition. Circulation 97, 1623-1631. doi:10.1161/01.CIR.97.16.1623

Verkerk, A. O. and Wilders, R. (2015). Pacemaker activity of the human sinoatrialnode: an update on the effects of mutations in HCN4 on the hyperpolarization-activated current. Int. J. Mol. Sci. 16, 3071-3094. doi:10.3390/ijms16023071

Verkerk, A. O., van Ginneken, A. C. G. andWilders, R. (2009). Pacemaker activityof the human sinoatrial node: role of the hyperpolarization-activated current, I(f ).Int. J. Cardiol. 132, 318-336. doi:10.1016/j.ijcard.2008.12.196

Viragh, S. and Challice, C. E. (1980). The development of the conduction system inthe mouse embryo heart. III. The development of sinus muscle and sinoatrialnode. Dev. Biol. 80, 28-45. doi:10.1016/0012-1606(80)90496-0

Wang, J., Vasaikar, S., Shi, Z., Greer, M. and Zhang, B. (2017). WebGestalt 2017:a more comprehensive, powerful, flexible and interactive gene set enrichmentanalysis toolkit. Nucleic Acids Res. 45, W130-W137. doi:10.1093/nar/gkx356

Waxman, J. S. and Yelon, D. (2009). Increased Hox activity mimics the teratogeniceffects of excess retinoic acid signaling. Dev. Dyn. 238, 1207-1213. doi:10.1002/dvdy.21951

Waxman, J. S., Keegan, B. R., Roberts, R. W., Poss, K. D. and Yelon, D. (2008).Hoxb5b acts downstream of retinoic acid signaling in the forelimb field to restrictheart field potential in zebrafish. Dev. Cell 15, 923-934. doi:10.1016/j.devcel.2008.09.009

Wen, Y. and Li, B. (2015). Morphology of mouse sinoatrial node and its expressionof NF-160 and HCN4. Int. J. Clin. Exp. Med. 8, 13383-13387.

Wiese, C., Grieskamp, T., Airik, R., Mommersteeg, M. T. M., Gardiwal, A., deGier-de, V. C., Schuster-Gossler, K., Moorman, A. F. M., Kispert, A. andChristoffels, V. M. (2009). Formation of the sinus node head and differentiation ofsinus nodemyocardium are independently regulated by tbx18 and tbx3.Circ. Res.104, 388-397. doi:10.1161/CIRCRESAHA.108.187062

Ye, W., Wang, J., Song, Y., Yu, D., Sun, C., Liu, C., Chen, F., Zhang, Y., Wang, F.,Harvey, R. P. et al. (2015). A common Shox2-Nkx2-5 antagonistic mechanismprimes the pacemaker cell fate in the pulmonary vein myocardium and sinoatrialnode. Development 142, 2521-2532. doi:10.1242/dev.120220

15

RESEARCH ARTICLE Development (2019) 146, dev173161. doi:10.1242/dev.173161

DEVELO

PM

ENT