Genome sequencing of adzuki bean (Vigna angularis)provides insight into high starch and low fataccumulation and domesticationKai Yanga,1, Zhixi Tianb,1, Chunhai Chenc,1, Longhai Luoc,1, Bo Zhaoa, Zhuo Wangc, Lili Yuc, Yisong Lia, Yudong Sunc,Weiyu Lia, Yan Chenc, Yongqiang Lia, Yueyang Zhangc, Danjiao Aia, Jinyang Zhaoc, Cheng Shanga, Yong Mac, Bin Wuc,Mingli Wanga, Li Gaoc, Dongjing Suna, Peng Zhangc, Fangfang Guoa, Weiwei Wanga, Yuan Lia, Jinlong Wanga,Rajeev K. Varshneyd,2, Jun Wangc,2, Hong-Qing Lingb,2, and Ping Wana,2

aCollege of Plant Science and Technology, Key Laboratory of New Technology in Agricultural Application, Beijing University of Agriculture, Beijing 102206,China; bState Key Laboratory of Plant Cell and Chromosome Engineering, Institute of Genetics and Developmental Biology, Chinese Academy ofSciences, Beijing 100101, China; cBeijing Genomics Institute-Shenzhen, Shenzhen 518083, China; and dInternational Crops Research Institute for theSemi-Arid Tropics, Hyderabad 502 324, India

Edited by Rod A. Wing, University of Arizona, Tucson, AZ, and accepted by the Editorial Board August 31, 2015 (received for review November 2, 2014)

Adzuki bean (Vigna angularis), an important legume crop, isgrown in more than 30 countries of the world. The seed of adzukibean, as an important source of starch, digestible protein, mineralelements, and vitamins, is widely used foods for at least a billionpeople. Here, we generated a high-quality draft genome sequenceof adzuki bean by whole-genome shotgun sequencing. The assem-bled contig sequences reached to 450 Mb (83% of the genome)with an N50 of 38 kb, and the total scaffold sequences were 466.7 Mbwith an N50 of 1.29 Mb. Of them, 372.9 Mb of scaffold sequenceswere assigned to the 11 chromosomes of adzuki bean by using asingle nucleotide polymorphism genetic map. A total of 34,183protein-coding genes were predicted. Functional analysis revealedthat significant differences in starch and fat content betweenadzuki bean and soybean were likely due to transcriptional abun-dance, rather than copy number variations, of the genes related tostarch and oil synthesis. We detected strong selection signals indomestication by the population analysis of 50 accessions includ-ing 11 wild, 11 semiwild, 17 landraces, and 11 improved varieties.In addition, the semiwild accessions were illuminated to have acloser relationship to the cultigen accessions than the wild type,suggesting that the semiwild adzuki bean might be a preliminarylandrace and play some roles in the adzuki bean domestication.The genome sequence of adzuki bean will facilitate the identificationof agronomically important genes and accelerate the improvement ofadzuki bean.

adzuki bean | genome sequence | starch synthesis genes | fat synthesisgenes | domestication

Adzuki bean (Vigna angularis var. angularis) was domesticatedin China ∼12,000 y ago (1) and is grown in more than 30

countries of the world, especially in East Asia (2, 3). Adzuki beanseed is an important source of protein, starch, mineral elements,and vitamins (4, 5). Because of its low caloric and fat content,digestible protein and abundant bioactive compounds, adzuki beanis referred to as the “weight loss bean” (6, 7). Given these char-acteristics, adzuki bean is widely used in a variety of foods (e.g.,paste in pastries, desserts, cake, porridge, adzuki rice, jelly, adzukimilk, ice cream) for at least a billion people (8). Furthermore,adzuki bean is a traditional medicine that has been used as adiuretic and antidote, and to alleviate symptoms of dropsy andberiberi in China (9, 10).Adzuki bean has broad adaptability, high tolerance to poor

soil fertility, and is a high-value rotation crop that contributes tothe improvement of soil condition by nitrogen fixation (11, 12).Additionally, adzuki bean can be used as a model species, espe-cially for non-oilseed legumes because of its short growth periodand small genome size (13, 14).

To accelerate the application of genomics for breeding and tobetter understand the biological mechanisms underlying its dis-tinct characters, we generated and analyzed the genome se-quence of adzuki bean variety Jingnong 6, which was bred by theBeijing University of Agriculture, China. The draft genome se-quence provides a resource for studying adzuki bean genomicsand genetics, promoting adzuki bean improvement.

Results and DiscussionGenome Sequencing and Assembly. The widely adopted Chinesecultivar “Jingnong 6” was targeted for whole genome shotgunsequencing by using the HiSeq 2000 sequencing platform. Intotal, 132.28 Gb sequence data were generated (SI Appendix,

Significance

Adzuki bean (Vigna angularis) is distinct in its high starch andlow fat accumulation. However, the underlying genetic basis isstill not well understood. In this study, we generated a high-quality draft genome sequence of adzuki bean by using whole-genome shotgun sequencing strategy. By comparative genomicand transcriptome analyses, we demonstrated that the signif-icant difference in starch and fat content between adzuki beanand soybean were caused by transcriptional abundance ratherthan copy number variations in the genes related to starch andoil synthesis. Furthermore, through resequencing of 49 pop-ulation accessions, we identified strong selection during do-mestication and suggested that the semiwild adzuki bean wasa preliminary landrace. Generally, our results provide insightinto evolution and metabolism of legumes.

Author contributions: Z.T., R.K.V., Jun Wang, H.-Q.L., and P.W. designed research; K.Y.,Z.T., B.Z., Yisong Li, Yongqiang Li, Y.Z., D.A., C.S., M.W., D.S., P.Z., F.G., W.W., Yuan Li,Jinlong Wang, H.-Q.L., and P.W. performed research; Jun Wang contributed new re-agents/analytic tools; C.H.C., L.H.L., Z.W., L.Y., Y.S., W.L., Y.C., J.Z., Y.M., B.W., and L.G.analyzed data; and K.Y., Z.T., C.H.C., R.K.V., H.-Q.L., and P.W. wrote the paper.

The authors declare no conflict of interest.

This article is a PNAS Direct Submission. R.A.W. is a guest editor invited by the EditorialBoard.

Freely available online through the PNAS open access option.

Data deposition: This Whole Genome Shotgun project has been deposited at DNA DataBank of Japan/European Molecular Biology Laboratory/GenBank under the accession no.JZJH00000000. The version described in this paper is version JZJH01000000. All short-readdata have been deposited into the Sequence Read Archive (SRA) under accession no.SRA259080. Raw sequence data of the transcriptomes have been deposited into theSRA under accession no. SRA260020.1K.Y., Z.T., C.H.C., and L.H.L. contributed equally to this work.2To whom correspondence may be addressed. Email: pingwan@bua.edu.cn, r.k.varshney@cgiar.org, wangj@genomics.org.cn, or hqling@genetics.ac.cn.

This article contains supporting information online at www.pnas.org/lookup/suppl/doi:10.1073/pnas.1420949112/-/DCSupplemental.

www.pnas.org/cgi/doi/10.1073/pnas.1420949112 PNAS | October 27, 2015 | vol. 112 | no. 43 | 13213–13218

AGRICU

LTURA

LSC

IENCE

S

Dow

nloa

ded

by g

uest

on

Mar

ch 2

3, 2

020

Table S1). After filtering out low-quality reads, 90.88 Gb high-quality sequences (168 × coverage of the adzuki bean genome)were obtained and used for assembly with SOAPdenovo software(15). Briefly, we generated 450 Mb contigs with an N50 size of38 kb, and 466.7 Mb scaffolds with an N50 size of 1.29 Mb (SIAppendix, Tables S2 and S3). The assembled scaffolds represent86.11% of adzuki bean genome (542 Mb) predicted by K-meranalysis (SI Appendix, Table S4). The GC content of adzuki beangenome is 34.8%, similar to other sequenced legume genomes(16–21) (SI Appendix, Fig. S1).To anchor the scaffolds to the chromosomes, we constructed a

high-density single-nucleotide polymorphism (SNP) genetic mapvia restriction-site associated DNA-sequencing (RAD-Seq) ap-proach by using 150 F2 individuals derived from a cross of cul-tivar Ass001 with the wild adzuki bean accession no. CWA108.The genetic map was comprised of 1,571 SNPs covering 11linkage groups, spanning 1,031.17 cM with an average of 4.33mapped SNPs per scaffold, and a mean marker distance of 0.67 cM(Fig. 1H and SI Appendix, Fig. S2). In total, 372.9 Mb of scaffolds,representing 79.9% of adzuki bean genome, were assigned to the11 pseudomolecules using these mapped SNPs (SI Appendix, Fig.S3 and Fig. 1B). A low conflict between genetic map and theassembled scaffolds indicates high quality of the genome assembly(SI Appendix, Table S5 and Fig. S4).

Repetitive Sequence Analysis and Gene Prediction. To analyze therepetitive sequences, we searched the genome sequence by usinga combination of de novo and homology tools, and found that∼44.51% (207.7 Mb) of adzuki bean genome consisted of re-petitive DNA. The percentage of the repetitive sequences inadzuki bean genome is higher than that of Medicago and Lotus,and lower than that of soybean, chickpea, and pigeonpea, butsimilar to common bean (SI Appendix, Table S6). Consistent withthe other legume genomes, the majority of the transposable elementswere retrotransposons (34.57% of genome), whereas DNA

transposons made up only 5.75% of the genome (SI Appendix,Table S7). Overall, ∼76.49% of the repetitive DNA was longterminal repeat retrotransposons of which 43.99% were Gypsy-type and 23.66% were Copia-type elements. An analysis of theadzuki bean genome with MIcroSAtellite (MISA) program (22)also allowed us to identify 16,230 simple sequence repeats(SSRs), containing di- (68.95%), tri- (18.12%), tetra- (1.98%),penta- (4.16%), and hexa- (6.79%) repeat units. Using theseinformation, we designed 9,038 SSR primer pairs, and found that24.7% of them exhibited polymorphism by analyzing 1,572 pairprimers (SI Appendix, Table S8). These SSRs will be a usefulresource for developing genetic markers and in genetic analysisof adzuki bean.To predict protein-coding genes, we sequenced the tran-

scriptomes of roots, stems, leaves, and seeds from three differentdevelopmental stages (beginning, full, and mature seed) of Jingnong6, and generated ∼309 million RNA reads, which were used formapping, clustering, and annotation of the whole genome (SI Ap-pendix, Table S9). The transcriptome assembly consisted of 59,909unitranscripts, of which 97.4% were covered by the genome as-sembly and 92.6% were captured in one scaffold with more than90% of the unitranscript length, further confirms the high quality ofthe genome assembly (SI Appendix, Table S10).Subsequently, we used a combination of de novo gene pre-

diction, homology-based search, and RNA-Seq to predict genemodels in the adzuki bean genome (SI Appendix, Fig. S5). A totalof 34,183 protein-coding genes were predicted (SI Appendix,Tables S2 and S11), of which 53.92% (18,432 genes) were sup-ported by RNA-Seq data, 81.33% (27,801 genes) by TrEMBL,59.67% (20,398 genes) by Swiss-Prot (23), and 63.35% (21,656genes) by Interpro (24) database (SI Appendix, Fig. S6). Therewere 6,196 genes (18.13% of total genes) that were functionallyunannotated (SI Appendix, Fig. S5). The genes were distributedunevenly with an increase in density toward the ends of thepseudomolecules (Fig. 1E). The density of expressed genes wassimilar to the pattern of gene distribution being lower or absentin the high repeat DNA regions of the pseudomolecules (Fig. 1D and F). Additionally, we identified 312 microRNAs, 307tRNAs, 3,730 rRNAs and 314 small nuclear RNAs (snRNA) (SIAppendix, Table S12).

Comparison of Gene Families with Other Sequenced Legume Genomes.Compared with six other sequenced legume genomes, the numberof predicted genes in the adzuki bean genome was lower thansoybean, pigeonpea, andMedicago, but higher than common bean,chickpea, and Lotus. Based on the proportion of the gene se-quence to the whole genome (total gene length/genome size),adzuki bean has a higher ratio (22.98%) than the other legumespecies with the exception of Medicago (SI Appendix, Table S11).For further analysis, 76,211 gene families of adzuki bean,

soybean, common bean, pigeonpea, and chickpea were clustered.We found that 12,582 gene families were common (orthologs) toall five legume genomes, whereas 827 gene families containing5,446 genes were specific to adzuki bean, more than soybean,common bean, and chickpea, and lower than that in pigeonpea(Fig. 2A). An analysis of these specific genes using Gene On-tology (GO) and Kyoto Encyclopedia of Genes and Genomes(KEGG) revealed that their functions mainly focused on thecategories of zinc ion binding, proteolysis, cysteine-type pepti-dase activity, two-component response regulator activity, andmitosis (Fig. 2B and SI Appendix, Table S13). We also found thatthe single copy gene orthologs in adzuki bean were significantlymore than those in soybean and similar to those in pigeonpea,chickpea, and common bean, whereas the multicopy orthologs inadzuki bean and the other three legume genomes were much lessthan that in soybean. This difference between adzuki bean andsoybean is most likely related to the additional whole genomeduplication (WGD) in soybean (17).In total, 3,508 adzuki bean transcription factor genes, be-

longing to 63 families, were identified by a sequence comparisonwith known transcription factors, and by searching for known

chr1

0 5 10 15

20

25

30

35

40

chr2

0

5

10

15

20

25

30

35

chr3

0

5

10

15

20

25

30

35

40

chr4

0

5

10

15

20

25

30

chr5

0

5

1015202530

chr6

051015

chr7

0

510

15

20

25

30

35

chr8

0

5

10

15

20

25

30

chr9

0

5

10

15

20

25

30 chr1

00

5

10

15

20

25

30

35

40chr11

0

5

10

15

20

25

85

0

45

0

45

0

35

0

1200

0HG

F

E

D

C

B

A

Fig. 1. Structure of the adzuki bean genome. (A) Gray lines connect duplicatedgenes among chromosomes of adzuki bean, and each line indicates a dupli-cation block containing at least 30 genes. (B) Pseudomolecules of adzuki bean.Distribution of GC content (C), repeat sequences density (D), gene density (E),expressive gene density (F), SSR density (G), and SNP density (H) in continuous200-kb windows were shown. The GC values range from 0.30 to 0.45.

13214 | www.pnas.org/cgi/doi/10.1073/pnas.1420949112 Yang et al.

Dow

nloa

ded

by g

uest

on

Mar

ch 2

3, 2

020

DNA-binding domains. These genes represented 10.26% of thetotal predicted genes in adzuki bean; this number was muchlower than those in common bean, soybean, and Lotus butsimilar to that in other sequenced legume species (SI Appendix,Table S14). In addition, we analyzed the genes encoding Rproteins in adzuki bean genome and compared them with otherlegume species (SI Appendix, Table S15). In total, 421 genescontaining NBS or LRR domains were detected in the adzuki beangenome, which is significantly lower than that in soybean andcommon bean. However, there were markedly more CC_NBSgenes in adzuki bean genome than that in soybean (SI Appendix,Table S15). This information should be helpful to identify the genesresponsible for plant disease and for disease resistance breeding.It has been reported that adzuki bean has quite low beany

flavor and contains digestible protein and abundant bioactivecompounds (4, 25, 26). Therefore, we analyzed genes related toflavonoid biosynthesis, lipoxygenase (LOX), and trypsin inhibitorin adzuki bean and other sequenced legume species, and nosignificant difference was observed in the gene number ratio (thegene number/the total gene number in the genome) (SI Appendix,Table S16). However, when we checked the expression of theLOX, which are responsible for the beany flavor in soybean (27),we found that their transcriptional amount was markedly lower inadzuki bean than that of soybean (SI Appendix, Fig. S7). Theseresults explain the low beany flavor of adzuki bean.Legume lectin, which is widely distributed among leguminous

plants, is a proteinaceous toxic factor that interacts with glyco-protein on the surface of red blood cells, causing them to ag-glutinate and a major antinutritional factor that inhibits theanimal growth, and affects the nutritional value and nutrientdigestibility (28). The lectin content (average 11.91 mg·g−1) ishigh in the grain of soybean (29), whereas it is low in adzuki bean(30). We found that legume lectin genes in adzuki bean showedsignificantly lower gene number ratio than that of other se-quenced legume species with an exception of chickpea. Consis-tently, the expression of the lectin genes [sum of the reads perkilobase of transcript per million reads mapped (RPKM) of in-dividual genes] was significantly lower in adzuki bean than thatof soybean, particularly for Le1, which is an important lectingene in soybean (gmx:100818710), but it is lacking in adzuki bean(SI Appendix, Fig. S8 and Table S17). This result indicates thatLe1 may play a main role in lectin accumulation in soybean seed.We analyzed gene family expansion and contraction (31) in

seven legumes and Arabidopsis. Among all 26,120 gene familiesof the eight species, 5.39% (1,407) and 7.83% (2,046) of the geneswere substantially expanded and contracted, respectively, in theadzuki bean during the 14 million years after speciation fromcommon bean, whereas soybean had many more expanded gene

families than common bean and adzuki bean (Fig. 3A and SIAppendix, Fig. S9), indicating that soybean has a universal ex-pansion in gene families, coincident with the larger gene number.We subsequently analyzed the gene functions in the expandedgene families of adzuki bean and found enrichment in GOcategories of zinc ion binding, proteolysis, cysteine-type pep-tidase activity, endopeptidase activity, endopeptidase inhibitoractivity, lipid binding, and lipid transport (SI Appendix, Tables S18and S19). For the contracted gene families of adzuki bean, theywere enriched in the categories of protein serine/threonine kinaseactivity, protein kinase activity, protein tyrosine kinase activity, anddefense response (SI Appendix, Table S20).

Genome Duplication and Synteny Analysis. Through ortholog search,we detected a total of 1,501 duplicated syntenic blocks within theadzuki bean genome. The gene number in a syntenic block rangedfrom 6 to 103 with an average of 11.7. The nucleotide diversity forthe fourfold degenerate third-codon transversion site (4DTv)showed one clear peak (4DTv∼0.36) in the adzuki bean genome(Fig. 3B), consistent with the whole genome duplication (WGD)event of Papilionoideae (32, 33). We did not identify the peak ofthe recent WGD (4DTv∼0.056) found in soybean, indicating thatadzuki bean, like most of the sequenced legumes, does not havethis glycine-specific event. Further phylogenetic analysis showedthat adzuki bean diverged from pigeonpea ∼19.0–32.5 millionyears ago (Mya), from soybean at ∼16.9–29.0 Mya, and fromcommon bean at ∼5.6–15.0 Mya (Fig. 3A).Synteny analysis revealed that, as expected, adzuki bean had

higher conservation with common bean than with other legumespecies (SI Appendix, Table S21 and Fig. 4A). Most of thechromosomes aligned between adzuki bean and common bean(e.g., adzuki bean chromosome 2 and common bean chromosome7, adzuki bean chromosome 5 and common bean chromosome 5,adzuki bean chromosome 1 and common bean chromosome 3, andadzuki bean chromosome 4 and common bean chromosome 9)(Fig. 4 A–C). However, some chromosomes of adzuki beanmatched to more than one chromosome of common bean, sug-gesting chromosome rearrangements had occurred in the two ge-nomes after speciation.We also did a comparison between adzuki bean and soybean.

These results showed that each adzuki bean chromosome matchedto several chromosomes of soybean, indicating that more ar-rangements occurred after speciation, and it is probably a result ofthe recent independent WGD in soybean (SI Appendix, Fig. S10).

Starch and Fatty Acid Biosynthesis and Metabolism Genes. The le-gume family, the second most important crop family, is divided intooil and non-oil types based on the storage compounds in the seed.Adzuki bean is a typical non-oil legume, whereas soybean belongs tothe oil type. Compared with soybean, adzuki bean seed has muchmore starch (57.06% vs. 25.3%) and much less crude fat (0.59% vs.

Pigeonpea

SoybeanAdzuki bean

Common bean

Chickpea

1448

780

354

827

21097

135

22

300

278

47

41

61254

42

147

362

12582

137

3041963

373

99

435

93

16

79134

108

490

57

A B

regulation of transcription, DNA -dependent

sequence-specific DNA binding transcription factor activity

sequence-specific DNA binding

zinc ion binding

proteolysis

cysteine-type pep�dase ac�vity

two -component response regulator activity

chromosome, centromeric region

mitosis

Fig. 2. Comparison of gene families in genomes of adzuki bean, commonbean, soybean, pigeonpea, and chickpea. (A) Venn diagram illustratingshared and unique gene families among the five species. (B) GO category ofunique gene families of adzuki bean.

A B

Fig. 3. Phylogeny of the seven genome-sequenced legumes. (A) Phyloge-netic tree and divergence time of genome-sequenced legume species, takingArabidoposis thaliana as an out group taxon. (B) Distributions of 4DTv dis-tance within adzuki bean, common bean, and soybean genomes.

Yang et al. PNAS | October 27, 2015 | vol. 112 | no. 43 | 13215

AGRICU

LTURA

LSC

IENCE

S

Dow

nloa

ded

by g

uest

on

Mar

ch 2

3, 2

020

22.5%) (5, 34). To investigate the bases behind this difference, weanalyzed genes related to starch and oil biosynthesis.Using the starch biosynthesis genes in rice (35, 36) as bait, we

performed an ortholog search in the adzuki bean and soybeangenomes. In total, fewer starch biosynthesis genes were found inadzuki bean than in soybean (27 vs. 46) (SI Appendix, Table S22),but the χ2 test result showed that the ratio of the starch bio-synthesis genes (number of starch synthesis genes/total genenumbers) was not significantly different between these two ge-nomes (P = 0.4135). We then examined the transcriptional ac-tivity of these genes at three seed development stages (beginning,full, and mature seed) of the two species with two biologicalreplicates collected in 2013 and 2014. The transcriptional level ofeach gene was normalized by using 47 constitutive expressiongenes in both species (SI Appendix, Table S23). We found thatboth the total transcriptional amount of starch biosynthesis genes(sum of the RPKM of individual genes) and average transcrip-tion of individual starch biosynthesis genes (mean of the RPKMof individual genes) in adzuki bean at the mature seed stage wassignificantly higher than that in soybean (Fig. 5A and SI Appendix,Figs. S11 and S12 and Table S24), although more starch bio-synthesis genes were detected in the soybean genome (SI Appendix,Table S22). However, no significant differences were observed atthe two earlier stages of seed development. Another interestingphenomenon was that the transcription of starch biosynthesis genescontinuously increased in adzuki bean, especially at the mature seedstage, whereas the expression of these genes in soybean decreasedat the full and mature seed stages (Fig. 5A and SI Appendix, Figs.S11 and S12).Subsequently, we searched for genes related to fatty acid

synthesis in plastids, the synthesis and storage of oil, and fattyacid degradation in the adzuki bean and soybean genomes. Al-though more genes were found in soybean than that in adzukibean, the copy number of the genes was not significantly dif-ferent relative to the total gene numbers of the two genomes(SI Appendix, Table S25). The transcription of the genes relatedto fatty acid synthesis in plastids and the synthesis and storage ofoil were markedly higher in soybean than those in adzuki bean(Fig. 5 B and C and SI Appendix, Figs. S11 and S12 and TablesS26 and S27). We also found that the two classes of genesexhibited differential expression patterns. Genes related to fattyacid synthesis in plastids exhibited a higher transcriptional level atearlier developmental stages in both adzuki bean and soybean,whereas the expression of oil synthesis and storage genes was

maintained in soybean but decreased in adzuki bean at the laterstage. The transcriptional amount of genes related to fatty aciddegradation in soybean was higher at the beginning seed stageand lower at the latter two stages than that observed in adzukibean (Fig. 5D and SI Appendix, Figs. S11 and S12 and TableS28). In addition, we checked the expression of five key genes thatencode the proteins that were involved in the conversion fromstarch to oil synthesis: amylase (an enzyme catalyzing the hydro-lysis of starch into sugars) (37), starch phosphorylase (an enzymeinvolving in phosphorolytic degradation of starch) (38), phos-phofructokinase (PFK, a kinase catalyzing the phosphorylationof fructose-6-phosphate to fructose-1,6-bisphosphate in the gly-colytic pathway) (39), glycerol-3-phosphate dehydrogenase (GPDH,an enzyme catalyzing the dihydroxyacetone phosphate to glycerol-3-phosphate, and serving as a major link between carbohydrate me-tabolism and lipid metabolism; ref. 40), and acetyl-CoA carboxylase(ACC, a biotin-dependent enzyme that catalyzes the irreversiblecarboxylation of acetyl-CoA to produce malonyl-CoA for the bio-synthesis of fatty acids) (41). All of the five genes showed a higherexpression in soybean at all three seed development stagescompared with adzuki bean (SI Appendix, Fig. S13 and TableS29). Based on these results, we speculated that the transcrip-tional amount of the genes related to starch and oil synthesisbetween soybean and adzuki bean are responsible for the dif-ferences of starch and fat content in the two species. Futuredetailed comparative genomic analysis of different species willhelp to answer the underlying mechanism.

Diversity and Domestication Analysis. In general, it is believed thatadzuki bean cultigen, including landraces and improved culti-vated variety, was domesticated from its putative wild progenitor(V. angularis var. nipponensis) (42). To investigate the geneticdiversity variation and selective sweeps during adzuki bean do-mestication, we performed whole genome resequencing of 11wild, 17 landraces, and 10 improved cultivars (SI Appendix, TableS30). Besides the wild type, landrace, and cultivar, another in-termediate type, semiwild adzuki bean (weedy), has been rec-ognized (2, 43). However, the genetic base of semiwild adzukibean was arguable. So far, it is not clear whether semiwild adzukibean was closely related to cultigen or wild adzuki bean or belongsto the landrace (44, 45). To make the resequencing populationrepresentative, we performed resequencing of 11 semiwild adzuki

Chr

1

0 1 2 3 4

Chr

2

0 1 23

Chr3

01

23

4 Chr40

12

3

Chr50

123

Chr601

Chr7

0123

Chr8

01

23

Chr9

01

23

Chr10

01

23

4Chr11

01

2

Chr11

012345

Chr

10

01234

Chr

09

0123

Chr08

01

23

4

5

Chr070

12

34

5Chr06 01

23

Chr0501234

Chr04

012

34

Chr03

01

23

4

5

Chr020

12

34

Chr01

0 1 2 3 4 5

Adzuki bean Chr2

Common bean Chr07

A B

C

Com

mon

bea

n

Adzuki bean

Common bean Adzuki bean

Chr1 Chr2 Chr3 Chr4 Chr5 Chr6Chr7 Chr8 Chr9

Chr10

Chr11

Chr1

Chr2

Chr3

Chr4

Chr5

Chr6

Chr7

Chr8

Chr9

Chr10

Chr11

Fig. 4. Genome synteny between adzuki bean and common bean chro-mosomes. Dotplot (A) and circos (B) figures of gene syntenic blocks (≥30 genesin each block) between adzuki bean and common bean chromosomes. Thechromosome orientations in B are clockwise. (C) Colinearity of genes of adzukibean chromosome 2 with common bean chromosome 7. The green and redlines indicate forward and reverse syntenic relationships, respectively.

A C D

B

Fig. 5. Comparison of the expression of the genes related to starch, fattyacid, and oil biosynthesis in adzuki bean and soybean. The transcriptionalamount of genes related to starch biosynthesis (A), synthesis of fatty acids inplastids (B), synthesis and storage of oil (C), and fatty acid degradation (D) inadzuki bean and soybean seeds at three developmental stages.

13216 | www.pnas.org/cgi/doi/10.1073/pnas.1420949112 Yang et al.

Dow

nloa

ded

by g

uest

on

Mar

ch 2

3, 2

020

beans as well. Therefore, a total of 49 adzuki bean accessions wereresequenced in this study.The whole genome resequencing yielded 254.41 Gb of high-

quality clean data. The sequencing depths of the 49 accessiongenomes ranged from 5.3× to 27.34×. The reads were aligned toadzuki bean reference genome by using SOAP2 (46) software,showing a mapping rate varied from 80.11 to 95.37%. In total,5,539,411 SNPs were identified among the 50 accessions acrossthe entire genome by SOAPsnp (47). Similar results were obtainedby BWA (48) along with SAMtools (49) analysis (SI Appendix,Table S31).The neighbor joining tree illustrated that the 50 adzuki bean

accessions could be divided into two main groups (Fig. 6A). The11 wild accessions clustered clearly together in a group, whereasthe remaining 39 accessions, including the semiwild accessions,landraces, and improved varieties, were grouped in anothergroup. The 11 semiwild accessions spread among the landracesand improved varieties. These results demonstrated that the semiwildadzuki bean had a closer relationship to cultigen adzuki beanthan that of wild type. Principal components analysis (50) also

suggested that the semiwild adzuki beans had a closer relation-ship with landraces and improved varieties than wild adzukibeans (Fig. 6B). The structure plots (51) showed the three typeswere separated into three groups when K was set as 3 withoutclear admixture (Fig. 6C). Based on all these comprehensiveanalyses, the semiwild adzuki bean seemed to be a distinct eco-type as the progenitor of cultigen accessions (44), rather than anescape from old cultivars and derivatives from a hybrid betweenwild and cultivars (42). Therefore, we classified the semiwildadzuki beans into landrace in the following selection analyses.We detected strong selection sweep signals [indicated by fix-

ation index (FST)] during the adzuki bean domestication, whichdistribute across all chromosomes (SI Appendix, Figs. S14 and S15and Table S31). The selection pressure detected between landracesand cultivars were obviously lower than that between wild and cul-tivars, which had much lower FST values, for example the FST valuesof the adzuki bean chromosome 1 (Fig. 6D). The results suggestedthat the domestication from the wild accession to cultivars is thor-ough and continuous. The genes in the selection regions enrichedmainly in the KEGG pathways of plant–pathogen interaction,plant hormone signal transduction, aminobenzoate degradation,cell cycle, and folate biosynthesis (SI Appendix, Fig. S16).

ConclusionIn this study, we obtained a high quality draft adzuki bean ge-nome sequence, in which more than 86% of the genome wasassembled, and approximately 80% of sequences were assignedto chromosomes. A total of 34,183 protein-coding genes werepredicted. Genome duplication analysis revealed that the adzukibean genome, unlike soybean, lacked the event of the recentwhole genome duplication. Compared with other sequenced legumegenomes, the adzuki bean genome had a higher synteny withcommon bean than that of soybean, pigeonpea,Medicago, chickpea,and Lotus. More interestingly, we revealed that the low fat and highstarch contents in adzuki bean seed were not caused by the genecopy number variation, but by gene expression amount in com-parison with soybean. Additionally, we also found that the semiwildadzuki bean would be a kind of preliminary landrace by populationanalysis with the 11 wild, 11 semiwild accessions, 17 landraces, and11 improved varieties, and detected strong selection signals in do-mestication. Generally, our results valuably reinforce the legumespecies genomes, provide insight into evolution and metabolicdifferences of legumes, and will accelerate studying of geneticsand genomics for adzuki bean improvement.

Materials and MethodsAdzuki bean (V. angularis var. angularis) cultivar Jingnong 6 (a popular va-riety in China) was used for whole genome sequencing. Forty-nine adzukibean accessions, including 11 wild adzuki bean (V. angularis var. nippo-nensis), 11 semiwild type, 17 landraces, and 10 improved cultivars wereperformed for their whole genome resequencing. The details of sequencing,assembly, annotation, genome analysis, and diversity and evolution analysiswere described in SI Appendix.

ACKNOWLEDGMENTS. We thank Drs. Huimin Wang (Beijing University ofAgriculture) and Yuncong Yao (Beijing University of Agriculture) for projectinitiation and support and Dr. Scott Allen Jackson (University of Georgia) forrevising the manuscript. Financial support was provided by Research Base andTechnological Innovation Platform Project of Beijing Municipal EducationCommittee Grant PXM2014_014207_000017, Enhancing Scientific ResearchLevel of the Beijing University of Agriculture Grant 1086716176, and Na-tional Natural Science Foundation of China Grants 31371694, 31272238, and31071474. This paper is dedicated to the late Professor Wenlin Jin, who de-voted his life to adzuki bean research. The sequencing Jingnong 6 variety wasbred by Prof. Jin.

1. Liu L, Bestel S, Shi J, Song Y, Chen X (2013) Paleolithic human exploitation of plant

foods during the last glacial maximum in North China. Proc Natl Acad Sci USA 110(14):

5380–5385.2. Tomooka N, et al. (2002) The Asian Vigna: Genus Vigna Subgenus Ceratotropis

Genetic Resources, ed Tomooka N (Kluwer Acad Pub, Dordrecht, The Netherlands),

pp 59–64.3. Kramer C, et al. (2012) Control of volunteer adzuki bean in soybean.Agri Sci 3(4):501–509.

4. Lin RF, et al. (2002)Minor Grain Crops in China, ed Yun L (Agric Sci Tech Press, Beijing),

pp 192–209.5. Xu N, et al. (2008) Establishment of an adzuki bean (Vigna angularis) core collection

based on geographical distribution and phenotypic data in China. Acta Agron Sin

34(8):1366–1373.6. Amarowicz R, et al. (2008) Antioxidant activity of extract of adzuki bean and its

fractions. J Food Lipids 15(1):119–136.

0 10 20 30 40

0.0

0.1

0.2

0.3

0.4

0.5

0.6

Chr1 (Mb)

Fst

Wild vs Landrace Wild vs Improved-cultivar Landrace vs Improved-cultivar

K=3

K=4

K=2

-0.6

-0.5

-0.4

-0.3

-0.2

-0.1

0

0.1

-0.1 0 0.1 0.2 0.9

Prin

cipal

com

pone

nt 2

Principal component 1

Improved-cultivarLandraceSemi-wild

Wild

0.2

1

Lxd3

ass001

AG110

Jingnong6Nongling5

KCA008KC

A013N

xd2

GM

652

Jxd3

Bao8

824-

17G

M90

4CW

A097

JWA09

2

CCA013

CCA031

CWA087

GM339

JCA009

CWA028

CWA085

JWA040JWA119JWA044KWA002KWA024KCA007

BD001

BD002

Suhong2

CCA047

CC

A050

CC

A053

AG86

CC

A052CW

A155

CCA0

65AG

57

CCA0

57

0.01

CWA098

CWA130

CWA143

CWA200CWA236

CWA221CWA194

CWA254

CWA246

NWA03CWA248

wil

d

cu

l t i ge

n

Landrace,Semi-wild

Improved-cultivar

A

B

D

C

Fig. 6. Genetic relationships in 50 adzuki bean accessions. (A) Neighbor joiningtree showing relationships of wild (blue), semiwild type (green), landraces(black), and improved varieties (red). (B) Principal component analysis of the 50adzuki bean individuals. (C) Structure analysis of the 50 adzuki bean individualsbased on SNP data. (D) Fst distribution on Chr1 for landrace vs. wild, im-proved variety vs. wild, and landrace vs. improved variety.

Yang et al. PNAS | October 27, 2015 | vol. 112 | no. 43 | 13217

AGRICU

LTURA

LSC

IENCE

S

Dow

nloa

ded

by g

uest

on

Mar

ch 2

3, 2

020

7. Kitano-Okada T, et al. (2012) Anti-obesity role of adzuki bean extract containingpolyphenols: in Vivo and in vitro effects. J Sci Food Agric 92(13):2644–2651.

8. Lestari P, et al. (2014) Genome-wide single nucleotide polymorphism discovery andvalidation in adzuki bean. Mol Breed 33(2):497–501.

9. Li S (1590, Ming Dynasty) (2010) Compendium of Materia Medica (Yunnan Educ Press,Kunming, China), p 255.

10. Chinese Pharmacopoeia Commission (2010) Chinese Pharmacopeia (the 2010 Edition), edsChinese Pharmacopoeia Commission (China Med Sci Press, Beijing), p 147.

11. Duan H (1989) Small bean. Edible Bean Crops, eds Duan H, Long J, Lin L, Xu H (Sci PubHouse, Beijing), pp 160–171.

12. Sikkema PH, et al. (2006) Response of adzuki bean to pre-emergence herbicides. Can JPlant Sci 86(2):601–604.

13. Parida A, Raina SN, Narayan RKJ (1990) Quantitative DNA variation between andwithin chromosome complements of Vigna species (Fabaceae). Genetica 82(2):125–133.

14. Yamada T, Teraishi M, Hattori K, Ishimoto M (2001) Transformation of azuki bean byAgrobacterium tumefaciens. Plant Cell Tissue Organ Cult 64(1):47–54.

15. Li R, et al. (2010) De novo assembly of human genomes with massively parallel shortread sequencing. Genome Res 20(2):265–272.

16. Schmutz J, et al. (2014) A reference genome for common bean and genome-wideanalysis of dual domestications. Nat Genet 46(7):707–713.

17. Schmutz J, et al. (2010) Genome sequence of the palaeopolyploid soybean. Nature463(7278):178–183.

18. Varshney RK, et al. (2012) Draft genome sequence of pigeonpea (Cajanus cajan), anorphan legume crop of resource-poor farmers. Nat Biotechnol 30(1):83–89.

19. Varshney RK, et al. (2013) Draft genome sequence of chickpea (Cicer arietinum)provides a resource for trait improvement. Nat Biotechnol 31(3):240–246.

20. Young ND, et al. (2011) The Medicago genome provides insight into the evolution ofrhizobial symbioses. Nature 480(7378):520–524.

21. Sato S, et al. (2008) Genome structure of the legume, Lotus japonicus. DNA Res 15(4):227–239.

22. Thiel T, Michalek W, Varshney RK, Graner A (2003) Exploiting EST databases for thedevelopment and characterization of gene-derived SSR-markers in barley (Hordeumvulgare L.). Theor Appl Genet 106(3):411–422.

23. Bairoch A, Apweiler R (2000) The SWISS-PROT protein sequence database and itssupplement TrEMBL in 2000. Nucleic Acids Res 28(1):45–48.

24. Zdobnov EM, Apweiler R (2001) InterProScan–an integration platform for the sig-nature-recognition methods in InterPro. Bioinformatics 17(9):847–848.

25. Yao Y, Cheng X, Wang L, Wang S, Ren G (2011) A determination of potential α-glucosidaseinhibitors from Azuki Beans (Vigna angularis). Int J Mol Sci 12(10):6445–6451.

26. Yamamoto A, et al. (1991) Partial purification and product specificity of adzuki beanlipoxygenase. J Jap Soc Food Sci Tech 38(3):214–216.

27. Davies CS, Nielsen SS, Nielsen NC (1987) Flavor improvement of soybean preparationsby genetic removal of lipoxygenase-2. J Am Oil Chem Soc 64(10):1428–1433.

28. Grant G (1989) Anti-nutritional effects of soyabean: A review. Prog Food Nutr Sci13(3-4):317–348.

29. Yang ML, Wang JA (2009) Agglutinin content determination and germplasm clusteranalysis in soybean seed. Soybean Bulletin 5:20–23.

30. Grant G, More LJ, McKenzie NH, Stewart JC, Pusztai A (1983) A survey of the nutri-tional and haemagglutination properties of legume seeds generally available in theUK. Br J Nutr 50(2):207–214.

31. De Bie T, Cristianini N, Demuth JP, Hahn MW (2006) CAFE: A computational tool for

the study of gene family evolution. Bioinformatics 22(10):1269–1271.32. Bertioli DJ, et al. (2009) An analysis of synteny of Arachis with Lotus and Medicago

sheds new light on the structure, stability and evolution of legume genomes. BMC

Genomics 10:45.33. Lavin M, Herendeen PS, Wojciechowski MF (2005) Evolutionary rates analysis of Le-

guminosae implicates a rapid diversification of lineages during the tertiary. Syst Biol

54(4):575–594.34. Jin WL, et al. (2006) Appraisal of quality characters for seeds of Chinese adzuki bean

local varieties. J Chin Cereals Oils Assoc 21(4):50–59.35. Tian Z, et al. (2009) Allelic diversities in rice starch biosynthesis lead to a diverse array

of rice eating and cooking qualities. Proc Natl Acad Sci USA 106(51):21760–21765.36. Ohdan T, et al. (2005) Expression profiling of genes involved in starch synthesis in sink

and source organs of rice. J Exp Bot 56(422):3229–3244.37. Manners DJ (1962) Enzymic synthesis and degradation of starch and glycogen. Adv

Carbohydr Chem 17:371–430.38. Tiwari R, Kumar A (2012) Starch phosphorylase: Biochemical and biotechnological

perspectives. Biotechnol Mol Bio Rev 7(3):69–83.39. Sweetlove LJ, Kruger NJ, Hill SA (2001) Starch synthesis in transgenic potato tubers

with increased 3-phosphoglyceric acid content as a consequence of increased

6-phosphofructokinase activity. Planta 213(3):478–482.40. Shen W, et al. (2010) Metabolic and transcriptional responses of glycerolipid path-

ways to a perturbation of glycerol 3-phosphate metabolism in Arabidopsis. J Biol

Chem 285(30):22957–22965.41. Abo Alrob O, Lopaschuk GD (2014) Role of CoA and acetyl-CoA in regulating cardiac

fatty acid and glucose oxidation. Biochem Soc Trans 42(4):1043–1051.42. Yamaguchi H (1992) Wild and weed azuki beans in Japan. Econ Bot 46(4):384–394.43. Yamaguchi H, Nikuma Y (1996) Biometric analysis on classification of weed, wild and

cultivated azuki beans. J Weed Sci Technol 41(1):55–62.44. Xu RQ, Tomooka N, Vaughan DA, Doi K (2000) The Vigna angularis complex: Genetic

variation and relationships revealed by RAPD analysis, and their implications for in

situ conservation and domestication. Genet Resour Crop Evol 47(2):123–134.45. Yoon MS, Lee J, Kim CY, Baek HJ (2007) Genetic relationships among cultivated and

wild Vigna angularis (Willd.) Ohwiet Ohashi and relatives from Korea based on AFLP

markers. Genet Resour Crop Evol 54(4):875–883.46. Li R, et al. (2009) SOAP2: An improved ultrafast tool for short read alignment.

Bioinformatics 25(15):1966–1967.47. Li R, et al. (2009) SNP detection for massively parallel whole-genome resequencing.

Genome Res 19(6):1124–1132.48. Li H, Durbin R (2009) Fast and accurate short read alignment with Burrows-Wheeler

transform. Bioinformatics 25(14):1754–1760.49. Li H, et al.; 1000 Genome Project Data Processing Subgroup (2009) The Sequence

Alignment/Map format and SAMtools. Bioinformatics 25(16):2078–2079.50. Price AL, et al. (2006) Principal components analysis corrects for stratification in ge-

nome-wide association studies. Nat Genet 38(8):904–909.51. Tang H, Peng J, Wang P, Risch NJ (2005) Estimation of individual admixture: Analytical

and study design considerations. Genet Epidemiol 28(4):289–301.

13218 | www.pnas.org/cgi/doi/10.1073/pnas.1420949112 Yang et al.

Dow

nloa

ded

by g

uest

on

Mar

ch 2

3, 2

020