• Bug reporting/known bugs
• File conversion/handling
• Coding style

The msatfinder script

1. Download the tar or zip file.
2. Unpack the downloaded file and cd to the new directory.
3. run

   ./msatfinder [options] file(s)

4. Er...
5. That's it.

You may make a file containing a list of all input files that you'd like msatfinder to run on, or supply their names on the command line (e.g. with a glob).

If you use a list, then run msatfinder as:

./msatfinder -l name_of_list_file

To run on the sample files provided you could provide their names or glob the file suffixes. For example:

./msatfinder *.gbk

Msatfinder can in fact be placed anywhere you like, for example in /usr/local/bin. As long as you have a configuration file in the same directory as your data then msatfinder will work. For example:

cd /home/user/msat_data/bacteria
/usr/local/bin/msatfinder *.gbk

Msatfinder should run without any need to configure it. If you'd like to change any of the parameters, please see the configuration section, below. There are a variety of options that can be used each time msatfinder is run, and these are described in the searching for microsatellites section.

Dependencies

• BioPerl.
• CGI.
• Config::Simple.
• File::Copy.
• Getopt::Std.
• Term::ReadLine.
• Term::ANSIColor.
• List::Util.

More information on installing Perl modules can be found at CPAN. However, many of them will be installed as standard with Perl 5.8.3 (if so, man or Perldoc should provide information).

External programs

The following external applications are used by msatfinder:

• EMBOSS.
• primer3. primer3_core should be installed on your PATH, as it is required for the eprimer3 part of EMBOSS to function (used by msatfinder to determine primers).

We use Gentoo Linux and BioLinux as the former has packages available for all these dependencies and the latter has them already installed.

Dependency search priorities

Msatfinder looks for (external program) dependencies as follows:

1. It will check whether the location given in the config file is present and executable.
2. If not, it will use “which” to look for the executable on one's path.
3. If none can be found, the script will die with a list of which dependencies were missing.

This system will allow you to specify alternative versions of dependencies, which will have priority over the ones on the user's PATH, if required. This may be useful if the user has more than one version of EMBOSS and running a specific version is required.

Overview

Msatfinder finds was designed to find perfect repeats (e.g. A(13) would be detected, but AAAAAATAAAAAA would not) in annotated (e.g. GenBank, EMBL, Swissprot) or unannotated (Fasta,raw) format files, but is also capable of finding interrupted microsatellites. It can be used to examine both protein and nucleic acid sequences. If given an annotated file it will extract information about each microsatellite and the sequence it is found in. In addition, for nucleic acid sequences it will determine whether it is possible to create PCR primers containing each repeat region found (EMBOSS/primer3 are required for this feature to work). The various features and output files may be controlled by editing the configuration file (msatfinder.rc) - notes are provided in the file on how to edit it, and more detail is given in this manual.

Input and output files

Input may be many separate sequence files, or multiple sequences in a file. To convert file formats readseq may be useful.

Input files

Input file types are automatically detected - if the file format cannot be determined a warning will be given. Allowed file types are:

• GenBank.
• Swissprot.
• EMBL.
• FASTA.
• ASCII/raw.

Input sequences may be amino acid or nucleic acid. Please note that if you use ‘ASCII’ format (i.e. raw sequence in a text file) then msatfinder will treat each new line as a new sequence, and blank lines are likely to cause it to fail.

You may name your input files as you please, but naming of the sequences in input files is particularly important for the use of msatfinder. As it is possible to use a file containing many sequences as the input to msatfinder, the script uses BioPerl to extract a unique name for each sequence. This is then used to name output files, such as the microsatellite FASTA files. It's therefore possible to have (for example) an input file called “MACO.gbk”, but the FASTA files will be called such things as “NC_004117.122038.CGA.6.fasta”. In the case of FASTA files, the unique identifier extracted is the first entry following the “>”. If you have no unique identifier in your files then msatfinder will attempt to generate one, but it is best if you can supply sequence files with each sequence already labeled.

Output files

When run, msatfinder creates some directories to store the output files — this is convenient as there are often large numbers of output files. If you install a version of msatfinder on your own computer then these directory names can be changed by editing the msatfinder.rc file. An html file (results.html) is created in the same directory in which msatfinder is run — open this in a browser to view all of msatfinder's output. Results.html contains links to the contents of the following seven subdirectories:

• Repeats: Contains various data files summarising genome and microsatellite information, e.g. taxonomy, number of microsatellites (see below). An index file summarising the results is created here.
• GFF: Contains GFF3 files describing the repeats. Some of these are equivalent to tab files, and can be viewed in Artemis. Click on the 'Artemis' link beside the file name. At the moment, this capability is only present when a single sequence has been supplied. NB. to view repeats, go to the bottom third of the Artemis window. Then scroll down this area. The repeats found by msatfinder are at the very end of the list (although their coordinates are scattered throughout the sequence). Please note that certain older genbank files (older than approx 1 year) will not work with the most recent Artemis release. To solve the problem, go to NCBI and download a newer version of your genbank file.
• Counts: A summary of all microsatellites found in a genome by length and type or exact motif. Unlike the index and count files in Repeats/, the files here preserve the exact motifs rather than just their content. For example TAT and AAT would be recorded as different microsatellites, whereas in some files they would be counted as the same as they have the same base content - this saves space. This directory also contains the thresholds suggested by Markov chain analysis or Poly, if these were chosen. If the former was used, another file containing the number of non-standard letter codes (nonstandard.count) is also generated.
• Msat_tabs, Flank_tabs: Contain feature tables for use with artemis. To use them, cd to either of these directories and run artemis file.gbk +file.tab. The msat_tab files show only the microsatellite itself, whilst the flank_tabs also include the flanking regions.
• Fasta: The sequence of the microsatellite plus flanking regions, in fasta format. These would be suitable for use in blast searches.
• MINE: These files contain the same information as shown in the “repeats” file in the Repeats directory (a summary of the details of each microsatellite), but one file per repeat is produced. These contain HTML formatting, and are designed to allow the creation of simple databases using our related MINE software, and are turned off by default.
• Primers: Primer files are produced for each microsatellite, if possible. By default it produces information for possible PCR primers. This option will be disabled if you enter a swissprot format amino acid sequence, but should be turned off if you intend to enter your AA sequence in FASTA format.

The following files will be found in the directory called “Repeats/”:

• sequence: This contains the information on each sequence, including number of microsatellites found, GC content, &c. A full listing of the column headers is shown here.
• repeats: This file contains the details of every individual microsatellite, plus similar genomic information to the sequence file. An explanation of the column headers is shown in the main msatfinder manual, here. Both this file and the sequence file may easily be imported into Open Office (or similar), or imported into a database.
• type.count: This is a table showing the number of microsatellites found, categorised by motif type (mono, di, tri &c.) and number of repeat units.
• motif.count: Similar to the type.count file, it shows the number of microsatellites found categorised by the bases/amino acids in the motif. These are ordered by the total base content only, thus AAT would be counted the same as TAT (exact summaries are available in the Counts/ directory).
• index: A handy summary of the results.
• primers.csv: A tabular summary of the information in all the primer files in the Primers/ directory.
• errors: file contains details of anything that looked unusual, e.g. very short sequences, features that did not match the sequence, &c.

Column headers for “repeats” and “sequence” files

The headers in these tables depend on the file format provided (e.g. Swissprot protein files don't have GC content). The headers are not found in the exact order shown below.

These column headings are found in both the repeats and sequence files in the Repeats/ directory.

sequence	The NC/Accession number/other identifier for this sequence.
specific_taxon, generic_taxon	Fields derived from the names of he directories you store your data in. Actual taxonomic information is captured in the fields below. You probably won't need this, but please e-mail if you need more information about this.
division,binomial,genus,species,strain/subspecies,common_name,organism,definition,taxid	Taxonomic features parsed from the genome file. These will not be available for some (e.g. FASTA) files.
strand	Whether the genome is single (ss) or double (ds) stranded.
alphabet	RNA, DNA or sometimes mRNA, depending on the genome annotation.
circular	Whether the genome is circular or linear
date	Submission date from the annotated file. These are reformatted into YYYY-MM-DD format.
specific_host,lab_host	This information on host range is sometimes provided in annotated files, but this field may often be blank.
notes	Any notes from the annotated file, if present.
genome_length	Length of the sequence.
no_of_coding_regions	The number of CDS annotations found in the genome.
total_nt_coding,percent_coding	The proportion of nucleotides that occur within CDS regions.

These additional column headings are found in the sequence file.

GC_content	GC content of the entire genome.
A, T, G, C	Base counts for the individual bases.
A/AT, C/GC	The ‘askew’ and ‘cskew’ values; the value is expressed as a percentage of no. of A in total number of AT, &c.
total_rep_length	Total length of all microsatellites in the genome.
pc_repeats	Microsatellites as a percentage of the genome.
msats, No_of_msats	The first value is set to 1 if any microsatellites were found, the second value is set to the total number found.

These additional column headings are found in the repeats file.

repeat	A unique name identifying this particular microsatellite.
flank_length	The length of the flanking regions (default: 300 each side).
repeat_plus_flank	Total length of the microsatellite and the flanks.
start, stop	Start and stop positions of the microsatellite within the genome.
motif_units	The motif and number of repeat units, e.g. AT(6).
motif, motifrevcom	The motif (e.g. AT) and its reverse complement (e.g TA).
repeat_units	Number of repeat units.
footprint	Total length of the microsatellite; no. of units multiplied by unit length.
dist_from_left, pc_from_left	Distance from the start of the sequence, as a number of bases and as a percentage, where the start of the microsatellite is found.
dist_from_right, pc_from_right	Distance from the end of the sequence, as a number of bases and as a percentage, where the start of the microsatellite is found.
motif_type	Mono, di, tri, tetra, penta or hexa.
reverse	whether the feature the msat is found in is on the forward or reverse strand (1 = reverse, 0 = forward or not in a feature).
primers	Was it possible to make a primer for this microsatellite (1 = yes, 0 = no).
GC_content(genome), GC_content(flank), GC_content(repeat)	GC content of the entire genome, the flanks only and the microsatellite itself.
coding_repeat	Set to 1 if msat is within a CDS region, 0 if it is not.
gene,product,protein	If the msat is withn a CDS region, msatfinder will parse all available annotations for that region and place them in these categories.

File names

Some of the output files produced by msatfinder have names such as “NC_000871.33400.AT.6.fasta”. These files are the fasta files of msat and flank sequence, primer files, mine files and feature tables.

The various sections are as follows.

NC_000871	NC number of the genome in which the microsatellite was found. This is derived from the unique identifier for each sequence (typically the accession number), extracted by msatfinder.
33400	Start position of the microsatellite.
AT	Repeat motif. This could be a combination of motifs (e.g. AT-TA) if the msat is interrupts
6	Number of repeat units.
fasta	Format of the file. Other suffixes include ‘primers.txt’ (primer files), ‘db’ (MINE files) and ‘count’ (summary of motif numbers). Feature tables have the filename NC_xxxxx.flank_tab or NC_xxxxx.msat_tab

Configuration

All of the parameters that can be customised by the user are found in the msatfinder.rc file. A brief description of how to set each of these is described in the file, but here we give some supplementary information. Variables in the COMMON and FINDER sections of the configuration file control the behaviour of msatfinder. Most of the default values will be acceptable in most cases.

• debug - if set to 1, will print extra debugging information. Set to 2 for even more.
• flank_size - the amount of sequence either side of the microsatellite that will be extracted and saved to the microsatellite FASTA file.
• mine_dir,repeat_dir,tab_dir,bigtab_dir,fasta_dir,prime_dir,align_dir,cont_dir - several variables that set the name of the subdirectories that will be created when the script is run.
• run_eprimer - set to 1 if you want to determine whether a primer can be made for each repeat.
• eprimer_args - the eprimer man page has more information on what to put here, if you are dissatisfied with the default (pick PCR primers). Please note that the “-task 0” option works with EMBOSS 2.8.0. If you have 2.9.0 then you should use “-primers” instead.
• eprimer - full path to the eprimer3 binary.
• primer3core - the full path to the primer3_core binary.
• override - turns off the following variables. It's easier than editing lots of lines in the config file.
  • artemis.
  • mine.
  • fastafile.
  • sumswitch.
  • screendump.
  • run_eprimer.
• motif_threshold - this is particularly important, as it defines the thresholds equal to or above which microsatellites will be detected, and which types will be detected. The types may be set to any length, and the lowest the thresholds can be set is 1, which will find every single base, pair of bases, triplet &c. It will take a long time to run if thresholds are set that low and the “regex” engine will not operate on such a small threshold. By default, mono-hexa will be searched for. Please refer to setting thresholds and motif types (below) for more information.
• artemis - turns on the Artemis feature tables.
• mine - turns on MINE summary files. These are equivalent to the "repeats" output file in the data they contain.
• fastafile - turns on whether or not a FASTA format file containing the sequence information for each microsatellite found will be generated.
• taxon information - two of the fields in the repeats and sequence files are “specific_taxon” and “generic_taxon”. See here.
• remote_link - used to put a hyperlink into MINE files for looking at the original genomes.
• sumswitch - determines whether or not the "repeats" output file will be created. This contains a large amount of information about each microsatellite and its genomic context, and can become rather large. However, it is very useful for importing into a database.
• screendump - prints out verbose information to the screen whilst running if set to 1.

Setting thresholds and motif types

It is possible to search for motif types of any length. However, configuring this may not be very intuitive. The default entries in msatfinder.rc dictate that msatfinder searches for motifs of length 1 to 6 (eg. A to ATTCCG). This, along with the thresholds for each motif length, is encoded by the line :

motif_threshold = "1,12|2,8|3,5|4,5|5,5|6,5"

In this case, a motif of length 1 (eg. A) must be repeated 12 times before being reported. A motif of length 2 (eg. AC) must be repeated 8 times, &c. If the user wishes to search for longer motifs, for example motifs 20 elements in length, then they simply need to add a pipe symbol (‘|’), followed by the motif length required and the corresponding threshold, eg. 2. The line in msatfinder.rc would then read:

motif_threshold = "1,12|2,8|3,5|4,5|5,5|6,5|20,2"

Likewise, if the user does not wish to search for a certain type of motif, eg. those of length 1, then they simply delete the '1,12' and the unneeded pipe symbol (pipes should not be at the beginning and end of the string). Doing so would reduce the time needed for searching.

Searching for microsatellites

Running msatfinder is described in the README and the installation section.

N.B.You must have an msatfinder.rc config file present in the same directory as your data for msatfinder to run (see the configuration section. If one is not present, the program will stop and give an error message suggesting that you place a copy of msatfinder.rc in the data directory. You may still use the help option if you don't have an msatfinder.rc file present (see below).

Msatfinder has a lot of output files, and users may sometimes want to back these up or suppress them. Running msatfinder with the --help or -h option will print out a list of the options available. The full list of options is:

Microsatellite searching engines

There are three “engines” currently implemented for microsatellite searching. These all operate in slightly different ways, and if you don't find the default useful then the others may prove effective. They are as follows.

• Regex (option 1): This uses fast regular expressions to search once through the sequence. It is the fastest method, but cannot detect very small repeats (threshold <3). This threshold is much lower than most people require, hence we have selected this engine as the default.
• Multipass (option 2): This steps through the sequence several times looking for microsatellites of one motif type on each pass, using regular expressions. This method may produce slightly different results, as it can find microsatellites that overlap.
• Iterative (option 3): This steps through the sequence one base at a time and attempts to construct a microsatellite at that position without using regular expressions. If it succeeds, it continues searching from the end of the last microsatellite. If you'd like to split the entire genome into repeat units this is the engine to use (it is slow, however).

The program is reasonably fast. For example, searching a large, microsatellite-rich bacterial genome such as Xylella fastidiosa (NC_002488, 2.68 Mb) with default thresholds takes under two minutes. The exact time will vary depending on your computing power.

Interrupted microsatellites

“Interrupted” microsatellites include several possible features. One is that a microsatellite tract consists of two motifs of the same class (e.g. dinucleotides) adjacent to each other (example 1). Another is that a long microsatellite tract may have one or more point mutations in it, making it appear to be several shorter tracts (example 2). Sometimes this latter category may include a pseudo-frameshift (example 3), thus appearing to be two different motifs.

Examples

1. ACACACACACACACAATGTGTGTGTGTGTGTG - a microsatellite tract consisting of (AC)8 with a single point mutation on the last unit (making AA), followed by a (TG)8 microsatellite.
2. AAAAAAAAAAGAAAAAAAAAA - what appears to be 2x (A)10 is in fact an (A)21, with a point mutation.
3. ATATATATATATTATATATATATAT - the insertion of an extra “T” into this msat makes it appear to be an AT followed by a TA motif, when it should be considered as an interrupted AT motif.

Once microsatellites have been detected by whichever engine has been selected, the complete list of msats is scanned to determine the relative positions of each microsatellite within each genome. Microsatellites are joined together when they meet two criteria:

• The distance from on microsatellite to the preceding one is equal to or less than the footprint of the current microsatellite.
• The current microsatellite and the preceding one are of the same motif length (mono-, di- &c.).

Each “cluster” of microsatellites thus found is combined into a single interrupted microsatellite, using the usual nomenclature. The motif type will be altered to contain all the motif types in the order they are found. So, example 1 (above) would become ac-tg.16, and example 2 would be a-a.21. In the latter case, the “-” is kept in to mark that this is an interrupted microsatellite.

The on-line version of msatfinder is available for those who don't want to do a local install. It offers the same features as the downloadable version, but due to server limitations can only accept sequences up to 70Mb. Explanations of each of the on-line options are shown below.

Motif selection

Threshold selection

ccccctccccctccccctccccct

Engines for finding perfect repeats

Advanced options

• sputnik - this uses the C program written by Chris Abajian (http://espressosoftware.com/pages/sputnik.jsp). It does not search for exact microsatellites, so that no motif is provided, and the microsatellite footprint is not equal to the motif length times the number of repeats.
• dust - used to filter nucleotide sequences for low-complexity areas.
• seg - finds areas of low compositional complexity in amino acid sequences, for example regions of biased amino acid composition like histidine-rich domains. Seg can be run with different options, including window (window size, default is 12), locut (low (trigger) complexity (default is 2.2)) and hicut (high (extension) complexity (default is 2.5)). More information for dust and seg available at http://www.molbiol.ox.ac.uk/analysis_tools/BLAST/BLAST_filtering.shtml#seg.

Download options

Find interrupted msats

If this box is checked (it's unchecked by default) then msatfinder will process the microsatellites found to determine if any could be joined into larger microsatellites, according to certain rules. Typically, about 10% of microsatellites found could be so joined.

Randomise sequence

Thresholds by Markov chain analysis

Thresholds by Poly

The third option for selecting microsatellite thresholds, along with manual choice and Markov chains, it to use the well-known program Poly, written by Jeff Bizzaro. The original publication is at http://www.biomedcentral.com/1471-2105/4/22. To speed up the calculations the algorithm has been transported from Python to C.

The user has two options:

1. Calculate Poly thresholds - Poly thresholds are calculated and stored in Fasta/sequenceID.poly.count, but not implemented
2. Implement Poly thresholds - Poly thresholds are calculated and implemented (ie. used to screen microsatellites)

Upload file

Choose a file from your local system to run msatfinder on. Files under 7Mb will be analysed and the results returned immediately to the user. Larger files require that the user submit their email address. Once the analysis is complete, a link to the results will be emailed to the user. Results are stored for 36 hours on the server before being deleted. If you have very large sequences, or many of them, we recommend a local installation and are happy to offer assistance, or run your sequence for you. See file formats for allowable file types.

Paste sequence

See file formats for input file formats. If you paste anything in here, it will be used instead of any uploaded file data, so make sure that this box is cleared if you'd like to upload a file.

Once you've submitted your sequence, click "search" and wait for a few moments (an animated picture will be shown whilst msatfinder runs). A brief summary will be displayed and a link to the downloadable file and the viewable results will be shown. These links will become inactive after a couple of hours, so please download them immediately if you'd like to keep the results.

Email

Output

The output will be viewable on-line for 36 hours, and can be downloaded as a tar or zip file within that time, so we recommend that you bookmark the link to your output. The various files that are included in the download directory are described below. The most important is "results.html" that includes links to all the other output files so that you may view them in a browser. Your input sequence and a configuration file will also be saved for future reference, in case of any problem with the results.

The output files produced by the online version of msatfinder are identical to those produced by the local version. A detailed description of these is available in the output file section. Please note that MINE files will not be produced by default - they must be enabled under “advanced options”.

What if it didn't work?

Occasionally, Msatfinder will fail to run on some input files. There are various reasons that this might happen, which are described below and also mentioned in the bugs section.

• The unique identifier extracted from the sequence may contain some characters that the operating system does not like. This does not happen very often, but if you suspect that this is the case then please send your sequence to Paul Swift. If you are using a FASTA format file you could also change the header to something inocuous like “> sequence1”.
• If your sequence contains lots of unknown bases then Msatfinder may have encountered a bug in Perl itself. The Perl developers know about this bug, but it seems that it cannot be fixed in the near future. Msatfinder will still work if the iterative engine is used, although this will be much slower.
• Another possibility is that your threshold settings are very low. As described here, low thresholds can cause ambiguities in the microsatellites found, and the program will not run if there is the potential for this. You can work around this problem by only looking for one type of microsatellite at a time, e.g. setting Msatfinder to find all monos of 3 or more units.

If you encounter a problem that does not seem to fit into these categories, please contact us and we will endeavour to fix the problem or analyse your data for you as soon as possible.

Bug reporting/known bugs

Though we are using this software successfully, there's a small chance of bugs turning up somewhere. Should you find any, please contact us (details below) and we will squash them. If you're using msatfinder, there are a few things that you ought to bear in mind.

• The number of interrupted microsatellites is dependent upon the thresholds chosen, as the program joins together microsatellites that are above the thresholds specified. An interrupted microsatellite whose individual parts are below the thresholds will not currently be detected.
• The program is reasonably fast, unless the thresholds are set very low or the genome is very large. The slowest part of the script is the section that parses genomic features, but this will not be needed if FASTA files are used.
• With low thresholds (e.g. 1), the counts may be skewed somewhat. For example, (AAT)7 could also be detected as seven (A)2 repeats. In general, it is not necessary to set thresholds this low, though. N.B. the default engine will not allow msats of less than three units to be detected, and will fail to run if ambiguities like this could be encountered, so this problem will only be seen if using one of the other engines.
• There is a problem in Perl's regex engine that can cause stack overflow in some cases. This generally means that if your data contains many unknown bases, msatfinder may die with a segmentation fault. Msatfinder is not the only application to be affected by this, and the Perl developers plan to change their code to overcome this problem. Meanwhile, if your data are affected, you can try a newer version of Perl or running msatfinder with the -e3 option (this is slow, but bypasses the problem).

We welcome suggestions for new features or other improvements.

File conversion/handling methods

For converting files one could use Seqret (from EMBOSS), or use readseq.

Coding style

The braces in this code are written thus:

This is simply because it's easier to read this way. However, if you disagree, you may wish to try this.