For the present study, the genomes were selected based on their quality and completeness. All genomes of the genus Streptomyces with status “Complete Genomes” were downloaded from Reference Sequence (Refseq; December 2019). After manual curation, 121 high quality genomes were included for the subsequent analysis. We further evaluated the genome assembly quality through the determination of genome completeness with BUSCO against the lineage dataset streptomycetales_odb10, which contains 145 species and 1,579 BUSCOs (Simão et al., 2015).

Genomes were downloaded from RefSeq in genbank and faa formats. We executed Roary to calculate the pan-genome for the genus Streptomyces (Page et al., 2015). Previously, we evaluated other software for pan-genomics studies such as BPGA (Chaudhari et al., 2016), GET_HOMOLOGUES (Contreras-Moreira and Vinuesa, 2013), Micropan in mode blast all vs all (Snipen and Liland, 2015) and Roary. The latter was selected because it had one of the lowest running times and generates the output for conservation of intergenic regions analysis, while producing similar results to the other tools. Roary requires genomes in gff3 format along with the sequences at the end of the file. To produce such files, we converted the genbank files using the BioPerl script bp_genbank2gff3.pl (McKay, 2004). The minimum percentage identity for BLASTp searches was set to 70%; the splitting of paralogs was blocked because it was required for the determination of the conserved IGRs. The maximum number of clusters was adjusted to 170,000. An alignment of core genes detected by Roary was created using MAFFT (Katoh and Standley, 2013) and we utilized this alignment to build a maximum likelihood phylogenetic tree using the software FastTree2 (Price et al., 2010) implemented in the Galaxy Europe server (Goecks et al., 2010).

In addition, we characterized the pan-genome of the genus Streptomyces based on the domain diversity of proteins encoded in the analyzed genomes. We used the R package Micropan version 1.2 (Snipen and Liland, 2015). Briefly, all the amino acid sequences of encoded proteins of the 121 genomes were annotated for their domain content with HMMER 3.3.1 (Eddy, 2011) against the Pfam-A database (Finn et al., 2014). Clustering made by Micropan is based on the presence of domains; thus, proteins sharing the same domains were grouped in the same gene family or cluster. The function BionomixEstimate of Micropan was implemented to extrapolate the size of the pan-genome using the presence/absence matrix resulting from both previous analyses. For both methodologies, we also determined the fluidity and the Jaccard distance for the genomes of the streptomycetes using the corresponding functions in Micropan.

All genes were classified as core, soft-core, shell, and cloud genes according to their presence among the genomes analyzed. Thus, genes present in the 121 strains were designated as core genes; genes present in more than 95% of strains (115 strains) were classified as soft-core genes; shell genes were those with a presence between 15 and 95% (19 and 114 strains), and genes present in less than 15% of the strains analyzed (less than 19 genomes) were assigned as cloud genes. For both methodologies, Roary and Micropan, we extracted representative sequences of the core, soft-core, shell, and cloud genes with in-house built Biopython scripts for subsequent functional annotation.

Genes resulting from Roary were translated into amino acid sequences; then, functional description of pan-genome categories defined for both methodologies was carried out determining gene ontology (GO) terms for the selected proteins (The Gene Ontology Consortium, 2019). This was performed with the Interproscan functional predictions of ORFs tool available in the Galaxy Europe Server (Quevillon et al., 2005). The results were summarized and plotted in WEGO 2.0 (Ye et al., 2018). Additional annotations were obtained through the WebMGA server for Clusters of Orthologous Groups (COG) assignments (Wu et al., 2011). Finally, the phylogenetic tree built with core genes along with information of the habitat and number of genes for each pan-genome category were visualized with Itol (Letunic and Bork, 2019).

In agreement with the phylogenetic tree, we defined three groups to analyze the conservation of IGRs in more closely related organisms; Streptomyces xiamenensis 318 and Streptomyces cattleya NRRL 8057 were left out of this analysis since no obvious relation with other Streptomyces was found. We estimated the conservation of intergenic regions across the streptomycetes using the software Piggy (Thorpe et al., 2018). The results of the previous analysis in Roary were used as input for Piggy. The software parameters nuc_id and len_id were set to 70, which is in accordance with the values used in Roary. Default values were used for the other parameters. Following this procedure, we analyzed the IGRs in the predefined groups of Streptomyces. The parameters used to analyze the groups of genomes were the same for the analysis of the complete set of genomes. Moreover, we aligned the IGRs conserved in more than 90% of the genomes included in each group against Rfam (version 14.5) database (Kalvari et al., 2021). Then, we explored for possible non-coding RNAs presence in these conserved IGRs with the software RNAz (Gruber et al., 2010). We previously filtered the IGRs alignments with the command rnazSelectSeqs.pl to preserve the sequences with a mean pairwise identity of 70%. Only the outputs with an overall RNA-class probability above 0.7 were considered as putative non-coding RNAs; their secondary structures were visualized with RNAfold (Lorenz et al., 2011) and their possible targets were defined using IntaRNA 2.0 (Mann et al., 2017).

The Galaxy wrapper of fastANI, with default parameters and using an all-versus-all genome comparisons, was implemented to calculate the average nucleotide identity (ANI) for the 121 selected strains (Jain et al., 2018). The heat map and dendrogram for the results of fastANI were generated using the libraries Seaborn and Matplotlib of Python (Hunter, 2007). The linkage method was the UPGMA algorithm and the pairwise distances between observations was the Euclidean metric. For those genomes with ANI values higher than 95% and with ambiguous taxonomic affiliations, we performed global genome alignments with progressiveMauve (Darling et al., 2010).

All 121 genomes were analyzed using ARTS 2.0 (available at https://arts.ziemertlab.com) with default settings. ARTS 2.0 used antiSMASH 5.1.1 for BGCs prediction (Blin et al., 2019). Since the lack of a proper methodology to define gene cluster boundaries, antiSMASH outputs a series of biosynthetic gene cluster regions; each region can be comprised by one or more co-localized “candidate” clusters; each “candidate” cluster defined by antiSMASH contains the biosynthetic machinery to produce a type of metabolite. In this work, we call BGC to each “candidate” cluster (for more information of antisSMASH definitions see: https://docs.antismash.secondarymetabolites.org/understanding_output/). After running the antiSMASH analysis, ARTS identifies BGCs co-localized with self-resistance enzymes (based on Resfam database), and with core genes (defined by ARTS using a database of actinomycetes genomes) with predicted horizontal gene transfer (HGT; Alanjary et al., 2017; Mungan et al., 2020). All predicted clusters of our interest were searched in the repository of the Minimum Information about a Biosynthetic Gene Cluster (MIBiG) database (Kautsar et al., 2020; available at https://mibig.secondarymetabolites.org).

The 3,750 genbank files, derived from the antiSMASH analysis, and the 33 genbank files corresponding to the BGCs prioritized by ARTS 2.0, were used as input for BGC similarity comparison using BiG-SCAPE 1.1.2 (Navarro-Muñoz et al., 2020). Analyses were made setting cutoff values at 0.3, 0.5, and 0.7 with and without the MIBiG database. Results of the networks, including the MIBiG database, were then filtered to remove comparisons between BGCs from the MIBiG database that did not display similarity with clusters from our analysis. Similarity comparison matrices were visualized using Cytoscape 3.8.2 (Shannon et al., 2003).

High quality genomes were included in the present investigation; apart from the status as “Complete genomes,” all genomes showed a high completeness and a reduced number of fragmented and missing genes (Supplementary Figure S1). Genome size ranges from 5.96 Mb for Streptomyces xiamenensis 318 to 12.01 Mb for S. hygroscopicus XM201; both strains also contain the minimum and the maximum protein coding genes with 5,100 and 9,385, respectively. The %G+C mean content is 71.77 +/− 0.81, which is an expected characteristic of members of the phylum Actinobacteria (Supplementary Figure S2; Dhakal et al., 2017). Most strains have a unique chromosome, although notably, the strain S. hygroscopicus limoneus KCTC 1717 has two chromosomes (Lee et al., 2016); all strains contain between one and four plasmids. Supplementary File S1 comprises all the metadata collected i.e., it contains the information of genome accession numbers, sequencing platform, coverage, and other genomic features such as the number of tRNAs and rRNAs in each genome.

We determined the pan-genome of the genus Streptomyces to establish their microbial diversity in terms of protein coding genes, domains content, and regulatory elements located in intergenic regions. For this purpose, we used different methodologies to accurately represent its entire gene repertoire. The analysis with Roary, to determine the diversity of protein coding genes, showed that the pan-genome of Streptomyces is clearly open (α<1, 0<γ<1) with a size of 145,462 clusters (Figures 1A,B). By using the BionomixEstimate function of Micropan, the current data allowed extrapolation to a total size of 273,372 clusters. These clusters were then classified according to their conservation level among the genomes analyzed. In concordance with this classification, we obtained 633 core genes, 1,080 soft-core genes, 6,040 shell genes, and 137,709 cloud genes; interestingly in the last group 81,568 were unique clusters, which means they were only present in one genome among all the considered strains.

Micropan estimated the pan-genome size of the genus Streptomyces as 8,973 protein families or clusters. Although, this number is clearly low compared to that determined by Roary, the power law fit, performed for the total size of the pan-genome and the number of new genes, showed that the pan-genome was still open (α<1, 0<γ<1), even though, it was set in the boundaries of a close pan-genome, as the value of the gamma parameter was close to zero (Figure 1C). The BionomixEstimate function, applied to both methods, displayed a similar core genome size of 600 and 589 for Roary and Micropan, respectively. Figure 1 summarizes the pan-genome calculations.

The fluidity of the pan-genome, which determines how dissimilar genomes are at a gene level, was estimated for both procedures employed to assess the genomic diversity of the genus Streptomyces (Kislyuk et al., 2011). The fluidity value was 0.53 +/− 0.099 for Roary and 0.22 +/− 0.031 for Micropan. This indicates that Streptomyces genomes differ 53%, on average, if the similarity of protein sequences are used to build the pan-genome, and 22% if their domain distributions are considered. A related assessment of genome diversity can be performed by the Jaccard distance distribution (Jaccard, 1912), which is roughly defined as one minus the number of genes shared by two genomes, divided by the total number of genes in these two genomes; the higher the value of Jaccard distance the more dissimilar the two genomes are. Overall, the Jaccard distance for both methodologies, Roary and Micropan, displayed similar distributions (Supplementary Figures S3A,B, respectively), centered at different mean values. Thus, highly similar genomes, as those genomes of the same species, possess the same genes/domains frequency giving values close to zero.

One of the main outcomes of a pan-genomic study is the determination of the genes shared by all members of a determined group, which corresponds to the core genome, previously defined. These core genes can be concatenated and aligned to define phylogenetic relationships among the members of a group, as this approach possesses higher resolution than using a single phylogenetic marker, e.g., 16S rRNA gene; thus, it has been suggested as the basis for bacterial phylogeny (Parks et al., 2018). Furthermore, a combination of results of multiple markers, such as the core genome phylogeny, above mentioned, and overall genome relatedness indices (from which ANI is the most broadly used), has been proposed to obtain precise taxonomic affiliations (Figueras et al., 2014). In this regard, we used both approaches, to explore the phylogeny in the genus Streptomyces.

Alignment of a set of 633 core genes, calculated by Roary, allowed the construction of a high confidence phylogenetic tree. The bootstrap values for all branches were above 0.9 being the majority equal to 1 (Figure 2). The number of genes in each genome, that are part of the different pan-genome categories, are also depicted in this figure, depending on the method used to determine the pan-genome. Based on the core genome phylogenetic tree, there was no clear relationship between the isolation source of the strains and its evolutionary relationship with other strains. Three clades were clearly distinguishable in the phylogenetic tree; they are highlighted in Supplementary Figure S4, for the sake of clarity.

Interestingly, S. hygroscopicus XM201 was set in group 1, while the other S. hygroscopicus strains were in group 3; in addition, the XM201 strain was closer to Streptomyces sp. 11-1-2 and S. violaceusniger Tu 4113, whose ANI values were 98.5 and 95.5, respectively. Other genomic features such as the genome size and the number of proteins encoded in the genome were also similar among these strains. Streptomyces lydicus strains formed a paraphyletic taxon; S. lydicus A02 was closer to S. gilvosporeus F607 than to the other strains classified as S. lydicus; this outcome was supported by ANI values in a range of 86–89%. A similar result was found for S. lydicus WYEC 108 whose ANI values were between 86 and 88% with other S. lydicus, and 96.7% with Streptomyces sp. NEAU-S7GS2. Lastly, Streptomyces sp. MOE7 contained an ANI value of 97.8% with S. lydicus GS93/23. ANI values of Streptomyces autolyticus CGMCC0516, Streptomyces malaysiensis DSM 4137 and Streptomyces sp. M56 were above 98% among them, which could indicate that they are the same species. Other strains that showed high ANI values between species and close relationship in the phylogenetic tree of core genes were: Streptomyces pratensis ATCC 33331 and Streptomyces sp. PAMC26508 (99.1 ANI); Streptomyces bacillaris ATCC 15855, Streptomyces sp. DUT11, Streptomyces sp. CFMR7 and Streptomyces sp. S8 (ANI values between 95.7 and 98.9%); both strains of Streptomyces globisporus with Streptomyces sp. 6063 and Streptomyces sp. Tue6075 (ANI values greater than 95.1); Streptomyces fradiae NKZ-259 and Streptomyces alfalfae ACCC40021 with an ANI value above of 99.9%, which might suggest they are the same strain, though further experimental studies are vital to prove it. The strains VK-A60T, KJ40, Fr-008, J1074, SM254, and SM17 all have ANI values greater than 95.8% among them.

An interesting clade is the one formed by Streptomyces sp. CB09001, the model organism Streptomyces coelicolor A3(2) and the biotechnologically important actinobacteria Streptomyces lividans TK24. A comparison of these strains showed that the genome of S. coelicolor A3(2) is almost 0.8 Mb larger than the S. lividans TK24 genome and 1.2 Mb larger than the one of Streptomyces sp. CB09001. Thus, we performed a global alignment of the genomes of these species to corroborate the observed relationship (Supplementary Figure S5). This assessment showed that the S. lividans TK24 genome is a reduced version of the S. coelicolor A3(2) genome, which has an additional region of about 0.6Mb in one of the telomeres. These results are in consensus with a recent study indicating that the ANI value between S. lividans and S. coelicolor is 99.0% (Vicente et al., 2018). Surprisingly, all the results observed in the core genome phylogenetic tree are validated by the corresponding ANI values among the species clustered together. In addition, a deep analysis for a possible reclassification of some species is suggested by the outcomes of the present study. These results can be observed in the Figure 2 and Supplementary Figure S4 for the phylogenetic tree and in the Supplementary Figure S6 for ANI values.

The number of IGRs was astonishingly high and variable, and no core group of IGRs was determined for the genus Streptomyces. The total number of IGR were 378,972; of these, 275,225 correspond to unique clusters of IGRs, which was more than twice the number of unique gene clusters obtained with Roary. As observed in Figure 3A, IGRs were only conserved across few strains. We further explored these results by analyzing the IGRs through the groups defined in the core genome phylogenetic tree, previously described (Supplementary Figure S4). In group 1, comprising 20 species (Figure 3B), the number of IGRs were still high (46,301 compared to 19,333 gene families), although, in this group there were 16 core IGRs. For group 2 (Figure 3C), which contained more species compared to group 1, only two IGRs were conserved in all 37 species belonging to this group; meanwhile, in group 3 (Figure 3D) integrated by 59 species, 131,213 clusters were unique IGRs and two were defined as core IGRs. Overall, IGRs clusters showed a pronounced drop as the number of genomes increased which differed from the behavior displayed for gene clusters, which were mainly unique or core genes (Figure 3A). Surprisingly, no annotations were retrieved from Rfam when representative sequences of these few conserved IGRs were searched. From these, 10, 1, 3 IGRs belonging to the groups previously defined, contain putative novel small or non-coding RNAs due to the conserved RNA secondary structures detected by RNAz (Supplementary File S2). The minimum free energy (MFE) structure of these predicted small RNAs (sRNAs) can be visualized in Supplementary Figures S7–S20. Additionally, we estimated the possible targets of these putative novel non-coding RNAs, in the genomes from which we extracted the representative IGRs sequences. Overall, we found six sequences that share full complementarity with the mRNA located down-stream, which suggests they can act as regulatory elements in the untranslated region of these genes; by other hand, multiple targets were detected that can interact with these sRNAs. The details of these analyses can be observed in the Supplementary File S2.

To investigate the IGRs conservation between more related Streptomyces species, we further analyzed the pan-genome of IGRs of S. lydicus (Figure 3E), Streptomyces clavuligerus (Figure 3F), Streptomyces albus (Figure 3G), and S. hygroscopicus (Figure 3H), as representatives of the three groups previously defined in the Streptomyces phylogeny. In the case of the paraphyletic group of S. lydicus, 248 IGR core clusters and 7,705 unique IGRs clusters were found. In S. clavuligerus, S. hygroscopicus, and S. albus, the number of core IGRs were 4,284, 4,597, and 4,706, respectively, which were considerably higher than the number of unique IGRs. This behavior agrees with the number of genes shared by the genomes, but it contrasts with the results obtained from the different groups of the phylogeny, and when all genomes of the genus were considered. Hence, we observed that IGRs are only conserved between phylogenetically related species.

Genes of the acquired pan-genome were then functionally classified. The COG functional enrichment demonstrated that the most conserved genes and family of proteins are those involved in primary metabolism and DNA processing functions (Figure 4). Interestingly, the abundance of secondary metabolism genes increases in less conserved genes, i.e., cloud genes. This tendency is more evident in the Micropan analysis, where protein domains of secondary metabolite genes represent more than 25% of total protein domains in cloud genes; lipid metabolism, frequently used for secondary metabolites production (Liu et al., 2013), also predominates.

The analysis of GO categories displayed similar results. Primary metabolism and catalytic processes as organic cyclic and heterocyclic compound binding are over-represented in core genes analyzed by Roary (Supplementary Figure S21). The GO enrichment in genes analyzed using Micropan evidences the abundance of genes involved in cellular and metabolic processes, as well as the abundance of the catalytic activity genes in all levels of conservancy, denoting the catalytic power of the genus.

Genomes were analyzed using ARTS 2.0 to prioritize BGCs more likely to produce an active metabolite, based on the presence of self-resistance enzymes co-localized within BGCs, as well as the presence of duplicated core genes with evidence of HGT (Alanjary et al., 2017; Mungan et al., 2020).

The analysis with antiSMASH displayed 3,750 regions of BGCs (Supplementary File S3). Since some BGCs can be co-localized in the same region (up to nine BGCs in one region), we separated the BGCs afterward to do the final count. However, it is worthy to point out that some co-localized BGCs could act as hybrid clusters, such as the modular system NRPS/T1PKS, which is widely found in the three domains of life (Wang et al., 2014).

Overall, 5,289 BGCs were identified in the 121 genomes analyzed, distributed in the 3,750 regions. Per order of frequency, non-ribosomal peptide synthetase (NRPS), terpene, type 1 polyketide synthase (T1PKS), and siderophore were the predominant BGC types, accounting for almost 50% of total predicted BGCs. Each genome accounts for 23–83 BGCs (average=44, median=42), Streptomyces griseochromogenes ATCC 14511 carried 83 BGCs and Streptomyces sp. CLI2509 carried 23. The biosynthetic potential of S. griseochromogenes ATCC 14511 was already unveiled by Wu et al. (2017a).

The set of Streptomyces strains analyzed carry 41 different types of BGCs out of 52 types defined by antiSMASH. The diversity in each genome goes between 10 and 26 types of BGCs (average and median=18). Streptomyces lydicus WYEC 108 was the strain that displayed the higher diversity, and Streptomyces koyangensis VK-A60T the lowest one. NRPS, terpene, and siderophore clusters were present in the 121 genomes (Figure 5; Supplementary Figure S4); although T1PKS and bacteriocin clusters were present in most of the strains, they were not found in Streptomyces exfoliatus A1013Y and Streptomyces xiamenensis 318, respectively. Furthermore, the ribosomally synthesized and post-translationally modified peptide (RiPP) clusters bottromycin and cyanobactin were only found in Streptomyces scabiei 87.22, and S. lydicus A02, respectively.

Although, no obvious relationship between the source of strains and their BGCs were found, there is a slight association between the frequency of BGC types and the genetic proximity (Supplementary Figure S4). For example, the clade of S. hygroscopicus displays similar frequency of NRPS, terpene, T1PKS, and siderophore; only the variety limoneus KCTC 1717 exhibited more bacteriocins in comparison to the varieties jinggangensis 5008 and the engineered jinggangensis TL01. In the case of S. lividans TK24 and S. coelicolor A3(2), both display a similar frequency of BGCs; yet, only S. lividans contains more terpenes in its genome. An interesting comparison is between Streptomyces sp. CNQ-509 and Streptomyces sp. WAC 06738; both strains are in the same clade but come from different isolation sources, marine and soil respectively, and mainly differ in the number of NRPS and T1PKS in their genomes.

The high BGCs variability in the genus was demonstrated with the cluster region comparison using BiG-SCAPE. This bioinformatic tool estimates distances between BGCs through the combination of the Jaccard index to determine the similarity of protein domains in the BGCs, the adjacency index that indicates the adjacent domains shared between BGCs, and the domain sequence similarity index, which calculate the sequence identity along with the domain copy number differences (Navarro-Muñoz et al., 2020). The network created by BiG-SCAPE using a cutoff of 0.3 – to identify interactions between BGCs producing similar compounds – displayed 2,359 nodes, and 12,969 edges (Figure 6A). A further comparison showed that 838 regions out of the 3,750 identified by antiSMASH are similar or have been already reported in the MIBiG database (Figure 6B; Supplementary File S4). Terpene, NRPS, siderophore, and ectoine are the clusters with the largest network similarity, whereas 1,204 cluster regions are unique within the analyzed genomes (Supplementary File S5).

To prioritize the search for antibiotics, ARTS uses BGC prediction from antiSMASH and displays the presence of self-resistance enzymes co-localized with BGC. In all 121 genomes analyzed, only 593 self-resistance genes were identified, distributed in 480 cluster regions out of the 3,750 regions identified by antiSMASH. On average, we identified five self-resistance genes in a genome; the maximum amount of self-resistance genes found in a genome was 12, in Streptomyces alfalfae ACCC40021. Streptomyces globisporus TFH56 was the unique strain without a self-resistance enzyme identified in its genome. Nevertheless, this strain can inhibit the growth of Botrytis cinerea, a gray mold pathogen that grows in tomato flowers (Cho and Kwak, 2019). Furthermore, we observed that NRPS and T1PKS are more frequently co-localized with self-resistance genes in comparison to other BGCs (Figure 5).

Another feature, considered in the prioritization of BGCs as possible producers of antibiotics, is the identification of core genes (defined by ARTS) within a biosynthetic cluster. Thus, 3,040 core genes were found in all genomes distributed in 1,490 regions. NRPS, terpene, and T1PKS offer the highest number of identified core genes (Figure 5). Additionally, core genes along with self-resistance genes were found in only 242 regions, being NRPS, T1PKS, and T2PKS the BGCs more frequently co-localized with both types of genes (Figure 5). Since some antibiotics target core genes, the producing bacteria tend to duplicate the gene and produce a homolog to avoid suicide. In this way, the presence of duplicated core genes in the BGC could lead to the prediction of the mode of action of the encoded antibiotic (Mungan et al., 2020). Applying a stricter filter to predict the antibiotic with its correspondent target, we only found 33 regions (distributed in 31 genomes) co-localized with self-resistance genes and core genes with evidence of duplication and HGT (Table 1); most of these clusters were NRPS (Figure 5).

In these 33 regions (Table 1), core genes were also classified as self-resistance genes; the diversity of these genes was low, presenting only three functions: glyceraldehyde-3-phosphate dehydrogenase (GAPDH) type I, proteasome, and DNA polymerase III β-subunit. Two self-resistance genes within the same region were found in only three genomes. In a NRPS-like cluster of Streptomyces lunaelactis MM109 the resistance targets were found in the C-terminal domain of GAPDH (GAPDH_C) and a metallo-β-lactamase, whereas in Streptomyces sp. WAC 01438 the T3PKS/NRPS/T2PKS cluster displayed two GAPDH_C as self-resistance genes, and Streptomyces sp. GGCR-6 presented a carboxyl transferase domain and GAPDH_C as resistance targets in a T1PKS cluster.

Using the approach of BGC prioritization, we identified BGCs with all elements needed to biosynthesize antibiotic molecules with a predicted mode of action. Some of the prioritized BGC display similarity with another prioritized cluster from a genetically related Streptomyces (Figure 6C), i.e., region 22 of S. lydicus 103 is similar to the region 16 of S. lydicus GS93 and the region 15 of Streptomyces sp. MOE7. Also, the region 5 of S. hygroscopicus XM201, the region 21 of Streptomyces violaceusniger Tu 4113 and the region 4 of Streptomyces sp. 11-1-2 are similar from each other. Likewise, the region 5 of S. autolyticus CGMCC0516 and the region 43 Streptomyces sp. M56 share sequence similarity. Intriguingly, the region 26 of Streptomyces collinus Tu 365 and the region 24 of S. cyaneogriseus noncyanogenus NMWT 1 are similar but both strains are not closely genetically related.

Of our prioritized BGCs only the region 19 of Streptomyces avermitilis MA-4680 is already described as the biosynthetic pathway of the antibiotic pentalenolactone, the region 18 also of S. avermitilis MA-4680 is reported in the MIBiG database (Kautsar et al., 2020) as a spore pigment cluster (although a bioactivity assay is not reported). Other four regions, along with the region 18 of S. avermitilis MA-4680, have similarity with a reported cluster in the MIBiG database (Table 1). Thus, we were able to perform a high throughput antibiotic screening using the bioinformatic tool ARTS 2.0, identifying interesting clusters that could be experimentally tested.