Bio::S...::msout(3pm)

💖 Bio::SeqIO::msout — the manual, but make it make sense 💖

okay bestie, gather round, because we're reading the manual for Bio::S...::msout(3pm) — User Contributed Perl Documentation, Bioperl edition — and honestly? it's giving niche-but-iconic. let's get into it. 🦋

---

✨ NAME ✨

Bio::SeqIO::msout — an input stream for output by Hudson's ms. yes that's the whole vibe. 💅

---

🎀 SYNOPSIS 🎀

do NOT come for this module directly, bestie. it is not for direct use. you use it through the Bio::SeqIO class, periodt.

---

🌸 DESCRIPTION 🌸

so ms is a whole thing — full citation: Hudson, R. R. (2002) Generating samples under a Wright-Fisher neutral model. Bioinformatics 18:337-8 — and you can find it living its best life at http://home.uchicago.edu/~rhudson1/source/mksamples.html. 🌈

r n, this object can read output from ms into seq objects. cute! but here's the tea ☕: bioperl does not have support for haplotypes created using an infinite sites model — that's the model where '1' flags a derived allele and '0' flags an ancestral allele, iykyk — so the sequences msout hands you back are coded using A, T, C and G instead.

to decode those bases, use the sequence conversion table (it's a hash, bestie) that you get back from get_base_conversion_table(). in that table, 4 and 5 get used when the ancestry is unclear. this should never ever happen when you're making files with ms, but it will get used when creating msOUT files from a collection of seq objects ( To be added later ). 🙈

alternatively — because we love options — use get_next_hap() to get a string with 1's and 0's instead of a seq object. slay. 🔥

🧚 Mapping to Finite Sites 🧚

not me telling you this object can now also map haplotypes created with an infinite sites model onto sequences of arbitrary finite length 👑. the details live in set_n_sites() — go read that for more. and a huge thank-you to Filipe G. Vieira <fgvieira@berkeley.edu> for the idea and the code. iconic behaviour, no cap. 💖

---

💌 FEEDBACK 💌

📬 Mailing Lists 📬

user feedback is an integral part of the evolution of this and other Bioperl modules — like, literally, it's the whole glow-up arc. so send your comments and suggestions, preferably to the Bioperl mailing list. your participation is much appreciated, fr fr. 🫶

bioperl-l@bioperl.org                  - General discussion
http://bioperl.org/wiki/Mailing_lists  - About the mailing lists

🐛 Reporting Bugs 🐛

found a bug? report it to the Bioperl bug tracking system so we can keep track of the bugs and their resolution — teamwork makes the dream work. bug reports can be submitted via the web:

https://github.com/bioperl/bioperl-live/issues

---

👑 AUTHOR - Warren Kretzschmar 👑

this module was written by Warren Kretzschmar. main character energy, honestly. 💜

email: wkretzsch@gmail.com

and plot twist — this module grew out of a parser written by Aida Andres. the lore is real. 🌟

---

⭐ COPYRIGHT ⭐

🌺 Public Domain Notice 🌺

this software/database is ``United States Government Work'' under the terms of the United States Copyright Act. it was written as part of the authors' official duties for the United States Government and thus cannot be copyrighted. this software/database is freely available to the public for use without a copyright notice. restrictions cannot be placed on its present or future use. it's giving public domain and we love that for it. 🎀

although all reasonable efforts have been taken to ensure the accuracy and reliability of the software and data, the National Human Genome Research Institute (NHGRI) and the U.S. Government does not and cannot warrant the performance or results that may be obtained by using this software or data. NHGRI and the U.S. Government disclaims all warranties as to performance, merchantability or fitness for any particular purpose. 😭

translation bestie: it's free, it's government work, and there is no warranty, periodt. don't come crying to NHGRI. 💀

---

✨ METHODS ✨

okay bestie, here's the part where we actually do stuff with the stream. every method, every arg, every warning — it's all here, no cap 💖

---

🧚 Methods for Internal Use

these are the little backstage gremlins. you can call them, but they're mostly doing their thing under the hood, periodt.

_initialize

• Title: _initialize
• Usage: $stream = Bio::SeqIO::msOUT->new($infile)
• Function: extracts basic information about the file.
• Returns: Bio::SeqIO object
• Args: no_og, gunzip, gzip, n_sites
• Details:
  • include the no_og flag if the last population of an msout file contains only one haplotype and you do NOT want that last haplotype to be treated as the outgroup (suggested when you're reading data created by ms). 👑
  • including n_sites (a positive integer) causes all output haplotypes to be mapped to a sequence of length n_sites. See set_n_sites() for the whole tea. 🫖

_read_start

• Title: _read_start
• Usage: $stream->_read_start()
• Function: reads from the filehandle $stream->{_filehandle} all information up to the first haplotype (sequence). Closes the filehandle if all lines have been read. 📖
• Returns: void
• Args: none

---

🌈 Methods to Access Data

this is the main character section. all the getters live here ✨

get_segsites

• Title: get_segsites
• Usage: $segsites = $stream->get_segsites()
• Function: returns the number of segsites in the msOUT file (according to the msOUT header line's -s option), or the current run's segsites if -s was not specified in the command line (in which case the number of segsites varies from run to run, lowkey chaotic).
• Returns: scalar
• Args: NONE

get_current_run_segsites

• Title: get_current_run_segsites
• Usage: $segsites = $stream->get_current_run_segsites()
• Function: returns the number of segsites in the run of the last read haplotype (sequence).
• Returns: scalar
• Args: NONE

get_n_sites

• Title: get_n_sites
• Usage: $n_sites = $stream->get_n_sites()
• Function: gets the number of total sites (variable or not) to be output.
• Returns: scalar if the n_sites option is defined at call time of new()
• Args: NONE
• Note:

WARNING: the final sequence length might not be equal to n_sites if n_sites is too close to the number of segregating sites in the msout file. 💀

set_n_sites

• Title: set_n_sites
• Usage: $n_sites = $stream->set_n_sites($value)
• Function: sets the number of total sites (variable or not) to be output.
• Returns: 1 on success; throws an error if $value is not a positive integer or undef
• Args: positive integer
• Note:

WARNING: the final sequence length might not be equal to n_sites if it is too close to the number of segregating sites. 😭

• n_sites needs to be at least as large as the number of segsites of the next haplotype returned
• n_sites may also be set to undef, in which case haplotypes are returned under the infinite sites model assumptions. 🦋

get_runs

• Title: get_runs
• Usage: $runs = $stream->get_runs()
• Function: returns the number of runs in the msOUT file (according to the msinfo line)
• Returns: scalar
• Args: NONE

get_Seeds

• Title: get_Seeds
• Usage: @seeds = $stream->get_Seeds()
• Function: returns an array of the seeds used in the creation of the msOUT file.
• Returns: array
• Args: NONE
• Details: in older versions, ms used three seeds. newer versions of ms seem to use only one (longer) seed. this function will return all the seeds found. ✨

get_Positions

• Title: get_Positions
• Usage: @positions = $stream->get_Positions()
• Function: returns an array of the names of each segsite of the run of the last read hap.
• Returns: array
• Args: NONE
• Details: the Positions may or may not vary from run to run depending on the options used with ms. 🙈

get_tot_run_haps

• Title: get_tot_run_haps
• Usage: $number_of_haps_per_run = $stream->get_tot_run_haps()
• Function: returns the number of haplotypes (sequences) in each run of the msOUT file ( according to the msinfo line ).
• Returns: scalar >= 0
• Args: NONE
• Details: this number should not vary from run to run. 👑

get_ms_info_line

• Title: get_ms_info_line
• Usage: $ms_info_line = $stream->get_ms_info_line()
• Function: returns the header line of the msOUT file.
• Returns: scalar
• Args: NONE

tot_haps

• Title: tot_haps
• Usage: $number_of_haplotypes_in_file = $stream->tot_haps()
• Function: returns the number of haplotypes (sequences) in the msOUT file. information gathered from the msOUT header line.
• Returns: scalar
• Args: NONE

get_Pops

• Title: get_Pops
• Usage: @pops = $stream->pops()
• Function: returns an array of population sizes (order taken from the -I flag in the msOUT header line). this array will include the last hap even if it looks like an outgroup. 💅
• Returns: array of scalars > 0
• Args: NONE

get_next_run_num

• Title: get_next_run_num
• Usage: $next_run_number = $stream->next_run_num()
• Function: returns the number of the ms run that the next haplotype (sequence) will be taken from (starting at 1). returns undef if the complete file has been read.
• Returns: scalar > 0 or undef
• Args: NONE

get_last_haps_run_num

• Title: get_last_haps_run_num
• Usage: $last_haps_run_number = $stream->get_last_haps_run_num()
• Function: returns the number of the ms run that the last haplotype (sequence) was taken from (starting at 1). returns undef if no hap has been read yet.
• Returns: scalar > 0 or undef
• Args: NONE

get_last_read_hap_num

• Title: get_last_read_hap_num
• Usage: $last_read_hap_num = $stream->get_last_read_hap_num()
• Function: returns the number (starting with 1) of the last haplotype read from the ms file
• Returns: scalar >= 0
• Args: NONE
• Details: 0 means that no haplotype has been read yet. it's reset to 0 every run. 🔄

outgroup

• Title: outgroup
• Usage: $outgroup = $stream->outgroup()
• Function: returns '1' if the msOUT stream has an outgroup. returns '0' otherwise.
• Returns: '1' or '0'
• Args: NONE
• Details: this method will return '1' only if the last population in the msOUT file contains only one haplotype. if the last population is not an outgroup then create the msOUT object using no_og as input flag. also, it returns 0 if the run has only one population. 🫶

get_next_haps_pop_num

• Title: get_next_haps_pop_num
• Usage: ($next_haps_pop_num, $num_haps_left_in_pop) = $stream->get_next_haps_pop_num()
• Function: first return value is the population number (starting with 1) the next hap will come from. the second return value is the number of haps left to read in the population from which the next hap will come.
• Returns: (scalar > 0, scalar > 0)
• Args: NONE

get_next_seq

• Title: get_next_seq
• Usage: $seq = $stream->get_next_seq()
• Function: reads and returns the next sequence (haplotype) in the stream
• Returns: Bio::Seq object or void if end of file
• Args: NONE
• Note: this function is included only to conform to convention. the returned Bio::Seq object holds a haplotype in coded form. use the hash returned by get_base_conversion_table() to convert 'A', 'T', 'C', 'G' back into 1, 2, 4 and 5. use get_next_hap() to retrieve the haplotype as a string of 1, 2, 4 and 5s instead. ✨

next_seq

• Title: next_seq
• Usage: $seq = $stream->next_seq()
• Function: alias to get_next_seq()
• Returns: Bio::Seq object or void if end of file
• Args: NONE
• Note: this function is only included for convention. it calls get_next_seq(). see get_next_seq() for details.

get_next_hap

• Title: get_next_hap
• Usage: $hap = $stream->next_hap()
• Function: reads and returns the next sequence (haplotype) in the stream. returns undef if all sequences in the stream have been read.
• Returns: haplotype string (e.g. '110110000101101045454000101'
• Args: NONE
• Note: use get_next_seq() if you want the haplotype returned as a Bio::Seq object. 👑

get_next_pop

• Title: get_next_pop
• Usage: @seqs = $stream->next_pop()
• Function: reads and returns all the remaining sequences (haplotypes) in the population of the next sequence. returns an empty list if no more haps remain to be read in the stream
• Returns: array of Bio::Seq objects
• Args: NONE

next_run

• Title: next_run
• Usage: @seqs = $stream->next_run()
• Function: reads and returns all the remaining sequences (haplotypes) in the ms run of the next sequence. returns an empty list if all haps have been read from the stream.
• Returns: array of Bio::Seq objects
• Args: NONE

---

💖 Methods to Retrieve Constants

base_conversion_table

• Title: get_base_conversion_table
• Usage: $table_hash_ref = $stream->get_base_conversion_table()
• Function: returns a reference to a hash. the keys of the hash are the letters 'A', 'T', 'G', 'C'. the values associated with each key are the value that each letter in the sequence of a seq object returned by a Bio::SeqIO::msout stream should be translated to.
• Returns: reference to a hash
• Args: NONE
• Synopsis:

# retrieve the Bio::Seq object's sequence
my $haplotype = $seq->seq;

# need to convert all letters to their corresponding numbers.
foreach my $base (keys %{$rh_base_conversion_table}){
        $haplotype =~ s/($base)/$rh_base_conversion_table->{$base}/g;
}

# $haplotype is now an ms style haplotype. (e.g. '100101101455')

---

and that's the whole page, bestie!! you made it all the way to the end, which is iconic of you honestly 🎀 this was perl v5.32.1, dated 2021-08-15, from the Bio::S...::msout(3pm) era of the perl docs — so like, a little vintage, but we stan a classic. go read some haplotypes for me, mwah 💅💜