Creating MAPMAKER/EXP data file

The MAPMAKER/EXP raw file, combined with the map file produced by OneMap, can be readily used for QTL mapping using R/qtl (Broman et al., 2008) or QTL Cartographer (Wang et al., 2010), among others.

Here, we briefly present how to set up this data file. For more detailed information see the MAPMAKER/EXP manual (Lincon et al., 1993), available here.

The first line of your data file should be:

data type xxxx

where xxxx is one of the following data types:

The second line should contain the number of individuals in the progeny, the number of markers, and the number of quantitative traits. So, for example, a valid line would be

10 5 2

for a data set with 10 individuals (yes, very small, but this is just an example), 5 markers, and 2 traits evaluated.

Then, the genotype information is included for each marker. The character * indicates the beginning of information of a marker, followed by the marker name. For instance, here is an example of such a file for an F2 population with 10 individuals, 5 markers, and 2 quantitative traits:

data type f2 intercross
10 5 2

*M1 A B H H A - B A A B
*M2 C - C C C - - C C A
*M3 D B D D - - B D D B
*M4 C C C - A C C A A C
*M5 C C C C C C C C C C

*weight 10.2 - 9.4 11.3 11.9 8.9 - 11.2 7.8 8.1 
*length 1.7 2.1 - 1.8 2.0 1.0 - 1.7 1.0 1.1

The codification for genotypes is the following:

The symbols option (not included in this example), used in MAPMAKER/EXP files, is also accepted (please, see its manual for details).

The quantitative trait data should come after the genotypic data and has a similar format, except the trait values for each individual must be separated by at least one space, a tab, or a line break. A dash (-) indicates missing data.

This file must be saved in plain text format using a simple text editor such as notepad on Microsoft Windows. Historically, MAPMAKER/EXP uses the .raw extension for this file; however, you can use any other extensions, such as .txt.

If you want to see more examples about this file type, open mapmaker_example_bc.raw and mapmaker_example_f2.raw, both available with OneMap and saved in the directory extdata on your computer, in the folder where you installed OneMap (use system.file(package="onemap") to see where it is located on your computer).

Now, let us load OneMap:

library(onemap)

To save your project anytime, type:

save.image("C:/.../yourfile.RData")

if you are using Windows, otherwise, adapt the code. Notice that you need to specify where to save and the name of the file. You can also use the toolbar, of course.

Creating OneMap data file

The OneMap data file has few differences compared to MAPMAKER/EXP format. As MAPMAKER/EXP format, the input OneMap file is a text file, where the first line indicates the cross-type and the second line provides information about the number of individuals, the number of markers, and the number of quantitative traits. Here, the format also supports keeping physical markers location information. The followed numbers indicate the presence/absence (1/0) of chromosome and position information and the presence/absence(1/0) of phenotypic data.

The third line contains sample IDs. Then, the genotype information is included separately for each marker. The character * indicates the beginning of information input for a new marker, followed by the marker name. Next, there is a code indicating the marker type according to:

Finally, after each marker type, comes the genotype data for the segregating population. Missing data are indicated with the character - (minus sign) and an empty space separates the information for each individual. Positions and phenotype information, if present, follows genotypic data with a similar structure. Details are found with the help of function read_onemap.

Here is an example of such file for 10 individuals, 5 markers (the two zeros in the second line indicate that there is no chromosome information, physical position information), and two phenotypic data, respectively). It is very similar to a MAPMAKER/EXP file, but has additional information about the cross_type.

data type f2 intercross
10 5 0 0 2
I1 I2 I3 I4 I5 I6 I7 I8 I9 I10
*M1 A.H.B  ab a - ab b a ab - ab b
*M2 A.H.B  a - ab ab - b a - a ab
*M3 C.A    c a a c c - a c a c
*M4 A.H.B  ab b - ab a b ab b - a
*M5 D.B    b b d - b d b b b d
*fen1 10.3 11.2 11.1 - 9.8 8.9 11.0 10.7 - 10.1
*fen2 42 49 - 45 51 42 28 32 38 40

The input file must be saved in text format, with extensions like .raw. It is a good idea to open the text files called onemap_example_f2.raw, onemap_example_bc, onemap_example_riself (available in extdata with OneMap and saved in the directory you installed it) to see how this file should be. You can see where OneMap is installed using the command system.file(package="onemap").

From MAPMAKER/EXP file

Once you created your data file with raw data, you can use OneMap function read_mapmaker to import it to OneMap:

mapmaker_example_f2 <- read_mapmaker(dir="C:/workingdirectory", 
                                     file="your_data_file.raw")

The first argument is the directory where the input file is located, modify it accordingly. The second one is the data file name.

In this example, an object named mapmaker_example_f2.raw was created. Notice that if you leave the argument dir blank, the file will be read from your current working directory. To set a working directory, see Introduction to R (Importing and Exporting Data).

mapmaker_example_f2 <- read_mapmaker(file= system.file("extdata/mapmaker_example_f2.raw", 
                                                       package = "onemap"))

For this example, we will use a simulated data set from an F2 population which is distributed along with OneMap. Because this particular data set is distributed along with the package, you can load it by typing:

data("mapmaker_example_f2")

To see what this data set is about, type:

mapmaker_example_f2

##   This is an object of class 'onemap'
##     Type of cross:      f2 
##     No. individuals:    200 
##     No. markers:        66 
##     CHROM information:  no 
##     POS information:    no 
##     Percent genotyped:  85 
## 
##     Segregation types:
##        AA : AB : BB -->  36
##         Not AA : AA -->  15
##         Not BB : BB -->  15
## 
##     No. traits:         1 
##     Missing trait values: 
##   Trait_1: 0

As you can see, the data consists of a sample of 200 individuals genotyped for 66 markers (36 co-dominant (AA, AB or BB), 15 dominant in one parent (Not AA or AA) and 15 dominant in the other parent (Not BB or BB) with 15% of missing data. You can also see that there is phenotypic information for one trait in the data set, that can be used for QTL mapping.

From OneMap raw file

The same procedure is made for the OneMap raw file, but, instead of using the function read_mapmaker we use read_onemap to read the OneMap format.

onemap_example_f2 <- read_onemap(dir="C:/workingdirectory", 
                                 inputfile = "your_data_file.raw")

In this example, an object named onemap_example_f2.raw was created. The data set containing the same markers and individuals of the mapmaker_example_f2.raw file. Would be a good idea to open these two files in a text editor and compare them to better understand the differences between the two kinds of input files. We can read the onemap_example_f2.raw using:

onemap_example_f2 <- read_onemap(inputfile= system.file("extdata/onemap_example_f2.raw", 
                                                        package = "onemap"))

Or, because this particular data are available together with the OneMap package:

data("onemap_example_f2")

To see what this data set is about, type:

onemap_example_f2

##   This is an object of class 'onemap'
##     Type of cross:      f2 
##     No. individuals:    200 
##     No. markers:        66 
##     CHROM information:  no 
##     POS information:    no 
##     Percent genotyped:  85 
## 
##     Segregation types:
##        AA : AB : BB -->  36
##         Not AA : AA -->  15
##         Not BB : BB -->  15
## 
##     No. traits:         1 
##     Missing trait values: 
##   Trait_1: 0

As you can see, the mean difference in the output object is that the read_onemap function keeps chromosome and position information. Because the objects mapmaker_example_f2 and onemap_example_f2 are practically the same, from now we will use only onemap_example_f2.

Importing data from the VCF file

If you are working with biallelic markers, as SNPs and indels (only codominant markers A.H.B), in VCF (Variant Call Format) files, you can import information from VCF to OneMap using onemap_read_vcfR function.

With the onemap_read_vcfR you can convert VCF file directly to onemap. The onemap_read_vcfR function keeps chromosome and position information for each marker at the end of the raw file.

We will use the same example file vcf_example_f2.vcf to show how it works.

Here we use the the vcfR package internally to help this conversion. The vcfR authors mentioned in their tutorials that RAM memory use is an important consideration when using the package. Depending of your dataset, the object created can be huge and occupy a lot of memory.

You can use onemap_read_vcfR function to convert the VCF file to onemap object. The parameters used are the vcf with the VCF file path, the identification of each parent (here, you must define only one sample for each parent) and the cross type.

vcf_example_f2 <- onemap_read_vcfR(vcf = system.file("extdata/vcf_example_f2.vcf.gz", package = "onemap"),
                                   parent1 = "P1", 
                                   parent2 = "P2", 
                                   cross = "f2 intercross")

Depending on your dataset, this function can take some time to run.

NOTE: From version 2.0.6 to 2.1.1005, OneMap had the vcf2raw function to convert vcf to .raw. Now, this function is defunct, but it can be replaced by a combination of onemap_read_vcfR and write_onemap_raw functions. See Exporting .raw file from onemap object session to further information about write_onemap_raw.

Graphical view of genotypes and allele depths (new!)

Function create_depth_profile generates dispersion graphics with x and y axis representing, respectively, the reference and alternative allele depths. The function is only available for biallelic markers in VCF files with allele counts information. Each dot represents a genotype for mks markers and inds individuals. If both arguments receive NULL, all markers and individuals are considered. Dots are colored according to the genotypes present in the OneMap object (GTfrom = onemap) or in the VCF file (GTfrom = vcf). An rds file is generated with the data in the graphic (rds.file). The alpha argument controls the transparency of the color of each dot. Control this parameter is a good idea when having a big amount of markers and individuals. The x_lim and y_lim control the axis scale limits, by default, it uses the maximum value of the counts.

Here is an example of the simulated dataset.

simu_f2_obj <- onemap_read_vcfR(vcf = system.file("extdata/vcf_example_f2.vcf.gz", package="onemap"), 
                                cross = "f2 intercross", 
                                parent1 = "P1", parent2 = "P2")
create_depths_profile(onemap.obj = simu_f2_obj,
                      vcfR.object =  system.file("extdata/vcf_example_f2.vcf.gz", package="onemap"), 
                      parent1 = "P1", 
                      parent2 = "P2", 
                      vcf.par = "AD", 
                      recovering = FALSE, 
                      mks = NULL, 
                      inds = NULL, 
                      GTfrom = "vcf", 
                      alpha = 0.1,
                      rds.file = "depths_f2.rds")

Selecting genotypes from vcf (GTfrom = “vcf”) the colors are separated by VCF genotypes: 0/0 homozygote for the reference allele, 0/1 heterozygote, and 1/1 homozygote for the alternative allele. Depending on your VCF, you can also have phased genotypes, which are presented by pipe (|) instead of the bar (/).

Using only recombinations information

mark_all_f2 <- make_seq(twopts_f2, "all")

class(mark_all_f2)

## [1] "sequence"

mrk_subset <- make_seq(twopts_f2, c(1, 3, 7))

mark_no_dist_f2 <- make_seq(twopts_f2, no_dist)

LGs_f2 <- group(mark_all_f2)

##    Selecting markers: 
##    group    1 
##     ............................................................
##     .............................

LGs_f2

##   This is an object of class 'group'
##   It was generated from the object "mark_all_f2"
## 
##   Criteria used to assign markers to groups:
##     LOD = 3 , Maximum recombination fraction = 0.5 
## 
##   No. markers:            90 
##   No. groups:             1 
##   No. linked markers:     90 
##   No. unlinked markers:   0 
## 
##   Printing groups:
##   Group 1 : 90 markers
##     M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 M12 M13 M14 M15 M16 M17 M18 M19 M20 M21 M22 M23 M24 M25 M26 M27 M28 M29 M30 M31 M32 M33 M34 M35 M36 M37 M38 M39 M40 M41 M42 M43 M44 M45 M46 M47 M48 M49 M50 M51 M52 M53 M54 M55 M56 M57 M58 M59 M60 M61 M62 M63 M64 M65 M66 SNP1 SNP2 SNP3 SNP5 SNP6 SNP7 SNP8 SNP9 SNP10 SNP11 SNP12 SNP13 SNP14 SNP15 SNP16 SNP17 SNP18 SNP19 SNP20 SNP21 SNP22 SNP23 SNP24 SNP26

(LGs_f2 <- group(mark_all_f2, LOD = LOD_sug, max.rf = 0.5))

##    Selecting markers: 
##    group    1 
##     ............................................................
##     .............................

##   This is an object of class 'group'
##   It was generated from the object "mark_all_f2"
## 
##   Criteria used to assign markers to groups:
##     LOD = 4.145858 , Maximum recombination fraction = 0.5 
## 
##   No. markers:            90 
##   No. groups:             1 
##   No. linked markers:     90 
##   No. unlinked markers:   0 
## 
##   Printing groups:
##   Group 1 : 90 markers
##     M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 M12 M13 M14 M15 M16 M17 M18 M19 M20 M21 M22 M23 M24 M25 M26 M27 M28 M29 M30 M31 M32 M33 M34 M35 M36 M37 M38 M39 M40 M41 M42 M43 M44 M45 M46 M47 M48 M49 M50 M51 M52 M53 M54 M55 M56 M57 M58 M59 M60 M61 M62 M63 M64 M65 M66 SNP1 SNP2 SNP3 SNP5 SNP6 SNP7 SNP8 SNP9 SNP10 SNP11 SNP12 SNP13 SNP14 SNP15 SNP16 SNP17 SNP18 SNP19 SNP20 SNP21 SNP22 SNP23 SNP24 SNP26

LGs_upgma <- group_upgma(mark_all_f2, expected.groups = 5, inter = F)
plot(LGs_upgma)

class(LGs_f2)

## [1] "group"

class(LGs_upgma)

## [1] "group.upgma"

set_map_fun(type = "haldane")

set_map_fun(type = "kosambi")

LG1_f2 <- make_seq(LGs_f2, 1)

LG1_f2

## 
## Number of markers: 90
## Too many markers  - not printing their names
## 
## Parameters not estimated.

class(LG1_f2)

## [1] "sequence"

LG1_rcd_f2 <- rcd(LG1_f2, hmm = FALSE)
LG1_rec_f2 <- record(LG1_f2, hmm = FALSE)
LG1_ug_f2 <- ug(LG1_f2, hmm = FALSE)

LG1_ser_f2 <- seriation(LG1_f2, hmm = F) # Will return an error (can not be used in this case)

LG1_mds_f2 <- mds_onemap(input.seq = LG1_f2, hmm = F)

## Using lodscore as value column: use value.var to override.
## Using lodscore as value column: use value.var to override.

## Stress: 0.179203004032195
## Mean Nearest Neighbour Fit: 47.5009683087036

ERROR: The linkage between markers 1 and 2 did not reach the OneMap default criteria. They are probably segregating independently

Linkage group 1

When possible (i.e., when groups have a small number of markers, in general up to 10 or 11), one should select the best order by comparing the multipoint likelihood of all possible orders between markers (exhaustive search). This procedure is implemented in the function compare. Although feasible for up to 10 or 11 markers, with 7 or more markers it will take a couple of hours until you see the results (depending of course on the computational resources available).

All our groups have more than 7 markers, so using the function compare is infeasible. Thus we will apply a heuristic that shows reliable results. First, we will choose a moderate number of markers, say 6, to create a framework using the function compare, and then we will position the remaining markers into this framework using the function try_seq. The way we choose these initial markers in inbred-based populations is somewhat different from what we did for outcrossing populations, where there is a mixture of segregation patterns (see the vignette for details).

In our scenario, we recommend two methods:

1. Randomly choose a number of markers and calculate the multipoint likelihood of all possible orders (using the function compare). If the LOD Score of the second-best order is greater than a given threshold, say, 3, then take the best order to proceed with the next step. If not, repeat the procedure.

2. Use some two-point based algorithm to construct a map; then, take equally spaced markers from this map. Then, create a framework of ordered markers using the function compare. Next, try to map the remaining markers, one at a time, beginning with co-dominant ones (most informative ones), then add the dominant ones.

You can do this procedure manually, in a similar way as done for outcrossing species (see the vignette for details). However, this procedure is automated in function order_seq, which we will use here (it will take some time to run):

LG1_f2_ord <- order_seq(input.seq = LG1_ug2_f2, n.init = 5,
                        subset.search = "twopt",
                        twopt.alg = "rcd", THRES = 3)

The first argument is an object of class sequence (LG1_ug2_f2). n.init = 5 means that five markers will be used in the compare step. The argument subset.search = "twopt" indicates that these five markers should be chosen by using a two-point method, which will be Rapid Chain Delineation, as indicated by the argument twopt.alg = "rcd". THRES = 3 indicates that the try_seq step will only add markers to the sequence which can be mapped with LOD Score greater than 3.

Check the order obtained by this procedure:

LG1_f2_ord # Results not shown in this vignette

Markers 5, 66, and 2 could not be safely mapped to a single position (LOD Score > THRES in absolute value). The output displays the safe order and the most likely positions for markers not mapped, where *** indicates the most likely position, and * corresponds to other plausible positions. (If you are familiar with MAPMAKER/EXP, you will recognize the representation).

To get the safe order, use:

LG1_f2_safe <- make_seq(LG1_f2_ord, "safe")

and to get the order with all markers (i.e., including the ones not mapped to a single position), use:

(LG1_f2_all <- make_seq(LG1_f2_ord, "force"))

which places markers 5, 66, and 2 into their most likely positions.

Although some old publications presented maps with only safe orders, we see no reason not to use the option force and recommend it for users. In the next, steps we can make modifications to the map based on more information.

The order_seq function can perform two rounds of the try_seq step, first using THRES and then THRES - 1 as the threshold. This generally results in safe orders with more markers mapped but takes longer to run. To do this, type (it will take some time to run):

LG1_f2_ord <- order_seq(input.seq = LG1_ug2_f2, n.init = 5,
                        subset.search = "twopt",
                        twopt.alg = "rcd", THRES = 3,
                        touchdown = TRUE)

## 
  |                                                                            
  |                                                                      |   0%
  |                                                                            
  |=                                                                     |   2%
  |                                                                            
  |==                                                                    |   3%
  |                                                                            
  |====                                                                  |   5%
  |                                                                            
  |=====                                                                 |   7%
  |                                                                            
  |======                                                                |   8%
  |                                                                            
  |=======                                                               |  10%
  |                                                                            
  |========                                                              |  12%
  |                                                                            
  |=========                                                             |  13%
  |                                                                            
  |==========                                                            |  15%
  |                                                                            
  |============                                                          |  17%
  |                                                                            
  |=============                                                         |  18%
  |                                                                            
  |==============                                                        |  20%
  |                                                                            
  |===============                                                       |  22%
  |                                                                            
  |================                                                      |  23%
  |                                                                            
  |==================                                                    |  25%
  |                                                                            
  |===================                                                   |  27%
  |                                                                            
  |====================                                                  |  28%
  |                                                                            
  |=====================                                                 |  30%
  |                                                                            
  |======================                                                |  32%
  |                                                                            
  |=======================                                               |  33%
  |                                                                            
  |========================                                              |  35%
  |                                                                            
  |==========================                                            |  37%
  |                                                                            
  |===========================                                           |  38%
  |                                                                            
  |============================                                          |  40%
  |                                                                            
  |=============================                                         |  42%
  |                                                                            
  |==============================                                        |  43%
  |                                                                            
  |================================                                      |  45%
  |                                                                            
  |=================================                                     |  47%
  |                                                                            
  |==================================                                    |  48%
  |                                                                            
  |===================================                                   |  50%
  |                                                                            
  |====================================                                  |  52%
  |                                                                            
  |=====================================                                 |  53%
  |                                                                            
  |======================================                                |  55%
  |                                                                            
  |========================================                              |  57%
  |                                                                            
  |=========================================                             |  58%
  |                                                                            
  |==========================================                            |  60%
  |                                                                            
  |===========================================                           |  62%
  |                                                                            
  |============================================                          |  63%
  |                                                                            
  |==============================================                        |  65%
  |                                                                            
  |===============================================                       |  67%
  |                                                                            
  |================================================                      |  68%
  |                                                                            
  |=================================================                     |  70%
  |                                                                            
  |==================================================                    |  72%
  |                                                                            
  |===================================================                   |  73%
  |                                                                            
  |====================================================                  |  75%
  |                                                                            
  |======================================================                |  77%
  |                                                                            
  |=======================================================               |  78%
  |                                                                            
  |========================================================              |  80%
  |                                                                            
  |=========================================================             |  82%
  |                                                                            
  |==========================================================            |  83%
  |                                                                            
  |============================================================          |  85%
  |                                                                            
  |=============================================================         |  87%
  |                                                                            
  |==============================================================        |  88%
  |                                                                            
  |===============================================================       |  90%
  |                                                                            
  |================================================================      |  92%
  |                                                                            
  |=================================================================     |  93%
  |                                                                            
  |==================================================================    |  95%
  |                                                                            
  |====================================================================  |  97%
  |                                                                            
  |===================================================================== |  98%
  |                                                                            
  |======================================================================| 100%

The output is too big to be included here, so please try it to see what happens. In short, for this particular sequence, the touchdown the step could additionally map markers 5 and 66, but this depends on the dataset. Let us continue our analysis using the order with all markers as suggested by the function order_seq:

(LG1_f2_final <- make_seq(LG1_f2_ord, "force"))

## 
## Printing map:
## 
## Markers           Position           Parent 1       Parent 2
## 
## 89 SNP24              0.00           a |  | b       a |  | b 
## 33 M33                5.53           a |  | b       a |  | b 
## 90 SNP26              6.97           a |  | b       a |  | b 
## 37 M37               13.33           a |  | b       a |  | b 
##  8 M8                19.15           a |  | o       o |  | o 
## 51 M51               30.87           o |  | a       o |  | o 
##  5 M5                34.62           o |  | a       o |  | o 
## 25 M25               37.24           a |  | b       a |  | b 
## 61 M61               42.79           a |  | b       a |  | b 
## 66 M66               44.47           a |  | b       a |  | b 
## 54 M54               51.42           a |  | b       a |  | b 
## 45 M45               53.93           a |  | o       o |  | o 
## 32 M32               58.95           o |  | o       o |  | a 
## 43 M43               63.28           a |  | b       a |  | b 
## 11 M11               66.80           a |  | o       o |  | o 
##  2 M2                69.58           a |  | o       o |  | o 
## 10 M10               77.02           a |  | b       a |  | b 
## 41 M41               85.50           a |  | b       a |  | b 
## 
## 18 markers            log-likelihood: -1256.273

rf_graph_table(LG1_f2_final)

Finally, to check for alternative orders, use the ripple_seq function:

ripple_seq(LG1_f2_final, ws = 5, LOD = 3)

The second argument, ws = 5, means that subsets (windows) of five markers will be permuted sequentially (5! orders for each window), to search for other plausible orders. The LOD argument means that only orders with LOD Score smaller than 3 will be printed.

The output shows sequences of five numbers, because ws = 5. They can be followed by an OK, if there are no alternative orders with LOD Scores smaller than LOD = 3 in absolute value, or by a list of alternative orders.

In this example, the first two windows showed alternative orders with LOD smaller than LOD = 3. However, the best order was that obtained with the order_seq function (LOD = 0.00). If there were an alternative order more likely than the original, one should check the difference between them and, if necessary, change the order with (for example) functions drop_marker (see Section about using an arbitrary order) and try_seq, or simply by typing the new order. For that, use LG2_f2_final$seq.num to obtain the original order; then make the necessary changes (by copying and pasting) and use the function map to reestimate the genetic map for the new order.

Linkage group 2

Here we will use the multipoint approach to re-order the LG2. The approach is the automatic usage of the try algorithm:

LG2_f2_ord <- order_seq(input.seq = LG2_ug_f2, n.init = 5,
                        subset.search = "twopt",
                        twopt.alg = "rcd", THRES = 3,
                        touchdown = TRUE, rm_unlinked = TRUE)

The second round of try_seq added markers 23, 29, 67, 44, 68, 36, 40, 26, 63, 31, 17, 12, 75, 74, 58, 35, 13, 6, 70, 72, 30, 69, 1, 46, 42, 3, 27 and 55 (try it; results not shown).

Get the order with all markers:

(LG2_f2_final <- make_seq(LG2_f2_ord, "force"))

## 
## Printing map:
## 
## Markers           Position           Parent 1       Parent 2
## 
## 23 M23                0.00           o |  | o       a |  | o 
## 60 M60                7.35           a |  | o       o |  | o 
## 29 M29               10.70           a |  | b       b |  | a 
## 67 SNP1              19.28           a |  | b       b |  | a 
## 44 M44               21.93           a |  | b       b |  | a 
## 68 SNP2              23.81           a |  | b       b |  | a 
## 36 M36               26.56           o |  | o       o |  | a 
## 40 M40               33.51           a |  | b       b |  | a 
## 26 M26               37.16           a |  | b       b |  | a 
## 63 M63               44.94           a |  | b       b |  | a 
## 34 M34               44.94           o |  | o       o |  | a 
## 31 M31               53.09           o |  | o       o |  | a 
## 17 M17               57.15           a |  | b       b |  | a 
## 73 SNP8              62.81           a |  | b       b |  | a 
## 12 M12               64.25           a |  | b       b |  | a 
## 75 SNP10             66.14           a |  | b       b |  | a 
## 74 SNP9              72.51           a |  | b       b |  | a 
## 58 M58               80.93           a |  | o       o |  | o 
## 35 M35               87.85           a |  | o       o |  | o 
## 13 M13               94.66           o |  | o       o |  | a 
##  6 M6                99.87           o |  | o       o |  | a 
## 71 SNP6             106.84           a |  | b       b |  | a 
## 70 SNP5             112.25           a |  | b       b |  | a 
##  7 M7               113.93           a |  | b       b |  | a 
## 72 SNP7             117.11           a |  | b       b |  | a 
## 30 M30              127.43           a |  | b       b |  | a 
## 69 SNP3             130.22           a |  | b       b |  | a 
##  1 M1               135.35           a |  | b       b |  | a 
## 46 M46              143.90           a |  | b       b |  | a 
## 53 M53              148.09           o |  | o       o |  | a 
## 42 M42              149.84           a |  | o       o |  | o 
##  4 M4               156.20           a |  | o       o |  | o 
##  3 M3               157.93           o |  | o       o |  | a 
##  9 M9               168.74           o |  | o       o |  | a 
## 27 M27              170.85           a |  | b       b |  | a 
## 55 M55              178.49           a |  | b       b |  | a 
## 
## 36 markers            log-likelihood: -2470.256

rf_graph_table(LG2_f2_final)

Our heatmap is already very good with an automatic approach, but I may want to see what happens if I remove a marker and try to reposition it. To remove a marker use function drop_marker:

LG2_edit <- drop_marker(LG2_f2_final, 23) # removing marker 23

After, you need to re-estimate the parameters using the same previous order. For that, use the map function:

LG2_edit_map <- map(LG2_edit) 

Warning: If you find an error message like:

Error in as_mapper(.f, ...) : argument ".f" is missing, with no default

It’s because the map function has a very common name and you can have in your environment other functions with the same name. In the case of the presented error, R is using the map function from purrr package instead of OneMap, to solve this, simply specify that you want the OneMap function with :: command from stringr package:

library(stringr)
LG2_edit_map <- onemap::map(LG2_edit) 

NOTE: (new!) If you are working with f2 intercross mapping population and have many markers (more than 60), we suggest to speed up map using BatchMap parallelization approach. See section Speed up analysis with parallelization for more information.

See what happened:

rf_graph_table(LG2_edit_map)

And try to reposition the marker:

(LG2_temp <- try_seq(input.seq = LG2_edit_map, mrk = 23))

## 
## LOD scores correspond to the best linkage phase combination
## for each position
## 
## The symbol "*" outside the box indicates that more than one
## linkage phase is possible for the corresponding position
## 
## 
##        Marker tested: 23
## 
##        Markers     LOD
##      =====================
##      |                   |
##      |             0.00  |  1  *
##      |  60               | 
##      |            -7.47  |  2  *
##      |  29               | 
##      |           -10.94  |  3  *
##      |  67               | 
##      |           -33.03  |  4  *
##      |  44               | 
##      |           -36.37  |  5  *
##      |  68               | 
##      |           -25.68  |  6  *
##      |  36               | 
##      |           -25.91  |  7  *
##      |  40               | 
##      |           -51.10  |  8  *
##      |  26               | 
##      |           -49.10  |  9  *
##      |  63               | 
##      |           -64.40  |  10  *
##      |  34               | 
##      |           -49.74  |  11  *
##      |  31               | 
##      |           -59.98  |  12  *
##      |  17               | 
##      |           -65.67  |  13  *
##      |  73               | 
##      |           -80.64  |  14  *
##      |  12               | 
##      |           -80.00  |  15  *
##      |  75               | 
##      |           -65.74  |  16  *
##      |  74               | 
##      |           -43.80  |  17  *
##      |  58               | 
##      |           -41.83  |  18  *
##      |  35               | 
##      |           -36.25  |  19  *
##      |  13               | 
##      |           -63.11  |  20  *
##      |   6               | 
##      |           -61.07  |  21  *
##      |  71               | 
##      |           -75.71  |  22  *
##      |  70               | 
##      |           -83.75  |  23  *
##      |   7               | 
##      |           -82.03  |  24  *
##      |  72               | 
##      |           -62.13  |  25  *
##      |  30               | 
##      |           -90.26  |  26  *
##      |  69               | 
##      |           -77.84  |  27  *
##      |   1               | 
##      |           -75.91  |  28  *
##      |  46               | 
##      |           -83.76  |  29  *
##      |  53               | 
##      |           -62.39  |  30  *
##      |  42               | 
##      |           -59.54  |  31  *
##      |   4               | 
##      |           -53.67  |  32  *
##      |   3               | 
##      |           -46.96  |  33  *
##      |   9               | 
##      |           -49.06  |  34  *
##      |  27               | 
##      |           -47.35  |  35  *
##      |  55               | 
##      |           -16.48  |  36  *
##      |                   |
##      =====================

The result shows us that the best position for this marker (with higher LOD) is the 1 (as before). Let’s now include the marker at this position using:

LG2_f2_final <- make_seq(LG2_temp, 1)

rf_graph_table(LG2_f2_final)

Check the final map (results not shown):

ripple_seq(LG2_f2_final, ws = 5)

Print it:

LG2_f2_final

## 
## Printing map:
## 
## Markers           Position           Parent 1       Parent 2
## 
## 23 M23                0.00           o |  | o       a |  | o 
## 60 M60               20.56           a |  | o       o |  | o 
## 29 M29               23.28           a |  | b       b |  | a 
## 67 SNP1              31.78           a |  | b       b |  | a 
## 44 M44               34.45           a |  | b       b |  | a 
## 68 SNP2              36.34           a |  | b       b |  | a 
## 36 M36               39.16           o |  | o       o |  | a 
## 40 M40               46.12           a |  | b       b |  | a 
## 26 M26               49.76           a |  | b       b |  | a 
## 63 M63               57.55           a |  | b       b |  | a 
## 34 M34               57.65           o |  | o       o |  | a 
## 31 M31               65.81           o |  | o       a |  | o 
## 17 M17               69.86           a |  | b       b |  | a 
## 73 SNP8              75.52           a |  | b       b |  | a 
## 12 M12               76.97           a |  | b       b |  | a 
## 75 SNP10             78.85           a |  | b       b |  | a 
## 74 SNP9              85.22           a |  | b       b |  | a 
## 58 M58               93.64           a |  | o       o |  | o 
## 35 M35              100.56           a |  | o       o |  | o 
## 13 M13              107.37           o |  | o       a |  | o 
##  6 M6               112.58           o |  | o       a |  | o 
## 71 SNP6             119.56           a |  | b       b |  | a 
## 70 SNP5             124.96           a |  | b       b |  | a 
##  7 M7               126.64           a |  | b       b |  | a 
## 72 SNP7             129.82           a |  | b       b |  | a 
## 30 M30              140.15           a |  | b       b |  | a 
## 69 SNP3             142.93           a |  | b       b |  | a 
##  1 M1               148.06           a |  | b       b |  | a 
## 46 M46              156.61           a |  | b       b |  | a 
## 53 M53              160.81           o |  | o       o |  | a 
## 42 M42              162.55           a |  | o       o |  | o 
##  4 M4               168.91           a |  | o       o |  | o 
##  3 M3               170.65           o |  | o       a |  | o 
##  9 M9               181.46           o |  | o       a |  | o 
## 27 M27              183.56           a |  | b       b |  | a 
## 55 M55              191.20           a |  | b       b |  | a 
## 
## 36 markers            log-likelihood: -2476.71

This is the final version of the map for this linkage group.

LG3_f2_ord <- order_seq(input.seq = LG3_ug_f2, n.init = 5,
                        subset.search = "twopt",
                        twopt.alg = "rcd", THRES = 3,
                        touchdown = TRUE)

(LG3_f2_final <- make_seq(LG3_f2_ord, "force"))

## 
## Printing map:
## 
## Markers           Position           Parent 1       Parent 2
## 
## 47 M47                0.00           a |  | b       a |  | b 
## 19 M19                7.68           a |  | o       o |  | o 
## 39 M39                8.98           o |  | o       o |  | a 
## 38 M38               15.59           a |  | b       a |  | b 
## 59 M59               23.86           a |  | o       o |  | o 
## 49 M49               24.32           o |  | o       a |  | o 
## 76 SNP11             32.50           a |  | b       a |  | b 
## 28 M28               33.81           a |  | b       a |  | b 
## 80 SNP15             36.11           a |  | b       a |  | b 
## 78 SNP13             42.82           a |  | b       a |  | b 
## 85 SNP20             52.38           a |  | b       a |  | b 
## 83 SNP18             57.63           a |  | b       a |  | b 
## 14 M14               59.25           a |  | b       a |  | b 
## 82 SNP17             61.56           a |  | b       a |  | b 
## 16 M16               67.72           a |  | b       a |  | b 
## 65 M65               72.93           a |  | b       a |  | b 
## 77 SNP12             78.35           a |  | b       a |  | b 
## 79 SNP14             88.64           a |  | b       a |  | b 
## 62 M62               98.09           a |  | b       a |  | b 
## 15 M15              101.64           a |  | o       o |  | o 
## 21 M21              105.74           a |  | b       a |  | b 
## 24 M24              107.72           o |  | o       a |  | o 
## 20 M20              114.56           o |  | o       a |  | o 
## 50 M50              138.46           o |  | o       a |  | o 
## 56 M56              141.91           a |  | b       a |  | b 
## 87 SNP22            153.60           a |  | b       a |  | b 
## 86 SNP21            159.91           a |  | b       a |  | b 
## 18 M18              163.36           a |  | b       a |  | b 
## 88 SNP23            164.83           a |  | b       a |  | b 
## 22 M22              170.26           a |  | b       a |  | b 
## 57 M57              173.68           a |  | b       a |  | b 
## 48 M48              179.82           a |  | b       a |  | b 
## 64 M64              185.51           a |  | o       o |  | o 
## 52 M52              187.45           a |  | b       a |  | b 
## 84 SNP19            189.98           a |  | b       a |  | b 
## 81 SNP16            196.26           a |  | b       a |  | b 
## 
## 36 markers            log-likelihood: -2633.912

rf_graph_table(LG3_f2_final)

LG3_edit <- drop_marker(LG3_f2_final, c(34,39,50,56, 20,24, 64))

## Warning in drop_marker(LG3_f2_final, c(34, 39, 50, 56, 20, 24, 64)): marker 34
## was not in the sequence

LG3_edit_map <- order_seq(LG3_edit) # We remove several markers maybe it's better to order again

## 
  |                                                                            
  |                                                                      |   0%
  |                                                                            
  |=                                                                     |   2%
  |                                                                            
  |==                                                                    |   3%
  |                                                                            
  |====                                                                  |   5%
  |                                                                            
  |=====                                                                 |   7%
  |                                                                            
  |======                                                                |   8%
  |                                                                            
  |=======                                                               |  10%
  |                                                                            
  |========                                                              |  12%
  |                                                                            
  |=========                                                             |  13%
  |                                                                            
  |==========                                                            |  15%
  |                                                                            
  |============                                                          |  17%
  |                                                                            
  |=============                                                         |  18%
  |                                                                            
  |==============                                                        |  20%
  |                                                                            
  |===============                                                       |  22%
  |                                                                            
  |================                                                      |  23%
  |                                                                            
  |==================                                                    |  25%
  |                                                                            
  |===================                                                   |  27%
  |                                                                            
  |====================                                                  |  28%
  |                                                                            
  |=====================                                                 |  30%
  |                                                                            
  |======================                                                |  32%
  |                                                                            
  |=======================                                               |  33%
  |                                                                            
  |========================                                              |  35%
  |                                                                            
  |==========================                                            |  37%
  |                                                                            
  |===========================                                           |  38%
  |                                                                            
  |============================                                          |  40%
  |                                                                            
  |=============================                                         |  42%
  |                                                                            
  |==============================                                        |  43%
  |                                                                            
  |================================                                      |  45%
  |                                                                            
  |=================================                                     |  47%
  |                                                                            
  |==================================                                    |  48%
  |                                                                            
  |===================================                                   |  50%
  |                                                                            
  |====================================                                  |  52%
  |                                                                            
  |=====================================                                 |  53%
  |                                                                            
  |======================================                                |  55%
  |                                                                            
  |========================================                              |  57%
  |                                                                            
  |=========================================                             |  58%
  |                                                                            
  |==========================================                            |  60%
  |                                                                            
  |===========================================                           |  62%
  |                                                                            
  |============================================                          |  63%
  |                                                                            
  |==============================================                        |  65%
  |                                                                            
  |===============================================                       |  67%
  |                                                                            
  |================================================                      |  68%
  |                                                                            
  |=================================================                     |  70%
  |                                                                            
  |==================================================                    |  72%
  |                                                                            
  |===================================================                   |  73%
  |                                                                            
  |====================================================                  |  75%
  |                                                                            
  |======================================================                |  77%
  |                                                                            
  |=======================================================               |  78%
  |                                                                            
  |========================================================              |  80%
  |                                                                            
  |=========================================================             |  82%
  |                                                                            
  |==========================================================            |  83%
  |                                                                            
  |============================================================          |  85%
  |                                                                            
  |=============================================================         |  87%
  |                                                                            
  |==============================================================        |  88%
  |                                                                            
  |===============================================================       |  90%
  |                                                                            
  |================================================================      |  92%
  |                                                                            
  |=================================================================     |  93%
  |                                                                            
  |==================================================================    |  95%
  |                                                                            
  |====================================================================  |  97%
  |                                                                            
  |===================================================================== |  98%
  |                                                                            
  |======================================================================| 100%

LG3_edit_map <- make_seq(LG3_edit_map, "force")

rf_graph_table(LG3_edit_map)

LG3_edit <- try_seq(LG3_edit_map, 34)
LG3_edit_temp <- make_seq(LG3_edit, 1) # Not included
LG3_edit <- try_seq(LG3_edit_map, 39)
LG3_edit_temp <- make_seq(LG3_edit, 3) 
LG3_edit_map <- LG3_edit_temp # include
LG3_edit <- try_seq(LG3_edit_map, 50)
LG3_edit_temp <- make_seq(LG3_edit, 22) # Not included
LG3_edit <- try_seq(LG3_edit_map, 56)
LG3_edit_temp <- make_seq(LG3_edit, 22) # Not included 
LG3_edit <- try_seq(LG3_edit_map, 20)
LG3_edit_temp <- make_seq(LG3_edit, 22) 
LG3_edit_map <- LG3_edit_temp # include
LG3_edit <- try_seq(LG3_edit_map, 24)
LG3_edit_temp <- make_seq(LG3_edit, 22) 
LG3_edit_map <- LG3_edit_temp # include
LG3_edit <- try_seq(LG3_edit_map, 64)
LG3_edit_temp <- make_seq(LG3_edit, 23) # Not included

LG3_f2_final <- LG3_edit_map

ripple_seq(LG3_f2_final, ws = 5)

idx <- which(LG3_f2_final$seq.num == 59) 
new_seq <- LG3_f2_final$seq.num
new_seq[idx:(idx+4)] <- c(59, 49, 76, 80, 28)
LG3_edit_seq <- make_seq(twopts_f2, new_seq)

LG3_edit_map <- onemap::map(LG3_edit_seq)

LG3_f2_final <- LG3_edit_map

rf_graph_table(LG3_f2_final)

Using the recombinations and the reference genome information

In our example, we have reference genome chromosome and position information for some of the markers, here we will exemplify one method of using this information to help build the genetic map.

With the CHROM information in the input file, you can identify markers belonging to some chromosome using the function make_seq with the rf_2pts object. For example, assign the string "1" for the second argument to get chromosome 1 makers. The output sequence will be automatically ordered by POS information.

CHR1 <- make_seq(twopts_f2, "1")
CHR1

## 
## Number of markers: 9
## Markers in the sequence:
## SNP1 SNP2 SNP3 SNP5 SNP6 SNP7 SNP8 SNP9 SNP10
## 
## Parameters not estimated.

CHR2 <- make_seq(twopts_f2, "2")
CHR3 <- make_seq(twopts_f2, "3")

CHR_mks <- group_seq(input.2pts = twopts_f2, seqs = "CHROM", unlink.mks = mark_all_f2,
                     repeated = FALSE)

##    Selecting markers: 
##    group    1 
##     ............................................................
##     ..............
##    Selecting markers: 
##    group    1 
##     ............................................................
##     ..................
##    Selecting markers: 
##    group    1 
##     ............................................................
##     .......

## There are one or more markers that grouped in more than one sequence

CHR_mks <- group_seq(input.2pts = twopts_f2, seqs = list(CHR1=CHR1, CHR2=CHR2, CHR3=CHR3), 
                     unlink.mks = mark_all_f2, repeated = FALSE)

CHR_mks

##   This is an object of class 'group_seq'
##   Criteria used to assign markers to groups:
##     LOD = , Maximum recombination fraction = 
## 
##   No. markers in input sequences:
##                        CHR1 :   9 markers
##                        CHR2 :   13 markers
##                        CHR3 :   2 markers
## 
##   No. unlinked input markers:   66 markers
## 
##   No. markers in output sequences:
##                        CHR1 :   9 markers
##                        CHR2 :   13 markers
##                        CHR3 :   2 markers
##   No. unlinked:                 0 markers
##   No. repeated:                 66 markers
## 
##   Printing output sequences:
##   Group CHR1 : 9 markers
##     SNP1 SNP2 SNP3 SNP5 SNP6 SNP7 SNP8 SNP9 SNP10 
## 
##   Group CHR2 : 13 markers
##     SNP11 SNP12 SNP13 SNP14 SNP15 SNP16 SNP17 SNP18 SNP19 SNP20 SNP21 SNP22 SNP23 
## 
##   Group CHR3 : 2 markers
##     SNP24 SNP26 
## 
##   Unlinked markers: 0 markers
##     
## 
##   Repeated markers: 66  markers
##      
##   Group CHR1 : 66 markers
##     M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 M12 M13 M14 M15 M16 M17 M18 M19 M20 M21 M22 M23 M24 M25 M26 M27 M28 M29 M30 M31 M32 M33 M34 M35 M36 M37 M38 M39 M40 M41 M42 M43 M44 M45 M46 M47 M48 M49 M50 M51 M52 M53 M54 M55 M56 M57 M58 M59 M60 M61 M62 M63 M64 M65 M66 
## 
##   Group CHR2 : 66 markers
##     M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 M12 M13 M14 M15 M16 M17 M18 M19 M20 M21 M22 M23 M24 M25 M26 M27 M28 M29 M30 M31 M32 M33 M34 M35 M36 M37 M38 M39 M40 M41 M42 M43 M44 M45 M46 M47 M48 M49 M50 M51 M52 M53 M54 M55 M56 M57 M58 M59 M60 M61 M62 M63 M64 M65 M66 
## 
##   Group CHR3 : 66 markers
##     M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 M12 M13 M14 M15 M16 M17 M18 M19 M20 M21 M22 M23 M24 M25 M26 M27 M28 M29 M30 M31 M32 M33 M34 M35 M36 M37 M38 M39 M40 M41 M42 M43 M44 M45 M46 M47 M48 M49 M50 M51 M52 M53 M54 M55 M56 M57 M58 M59 M60 M61 M62 M63 M64 M65 M66

CHR_mks$repeated

## $CHR1
##  [1]  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
## [26] 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50
## [51] 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
## 
## $CHR2
##  [1]  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
## [26] 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50
## [51] 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
## 
## $CHR3
##  [1]  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
## [26] 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50
## [51] 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

CHR_mks$sequences$CHR1

## 
## Number of markers: 9
## Markers in the sequence:
## SNP1 SNP2 SNP3 SNP5 SNP6 SNP7 SNP8 SNP9 SNP10
## 
## Parameters not estimated.

# or
CHR_mks$sequences[[1]]

## 
## Number of markers: 9
## Markers in the sequence:
## SNP1 SNP2 SNP3 SNP5 SNP6 SNP7 SNP8 SNP9 SNP10
## 
## Parameters not estimated.

Draw_map

We can draw a genetic map for all linkage groups using the function draw_map. First, we have to create a list of ordered linkage groups:

maps_list <- list(LG1_f2_final, LG2_f2_final, LG3_f2_final)

Then use function draw_map for this list:

draw_map(maps_list, names = TRUE, grid = TRUE, cex.mrk = 0.7)

We also can draw a map for a specific linkage group:

draw_map(LG1_f2_final, names = TRUE, grid = TRUE, cex.mrk = 0.7)

Function draw_map draws a very simple graphic representation of the genetic map. More recently, we developed a new version called draw_map2 that draws a more sophisticated figure. Furthermore, once the distances and the linkage phases are estimated, other map figures can be drawn by the user with any appropriate software.

Draw_map2

The same figures did with draw_map can be done with the draw_map2 function. But it has different capacities and arguments. Here are some examples, but you can find more options on the help page ?write_map2.

Let’s draw all three groups built:

draw_map2(LG1_f2_final, LG2_f2_final, LG3_f2_final, main = "Only linkage information", 
          group.names = c("LG1", "LG2", "LG3"))

You can include all sequence objects referring to the groups as the first arguments. The main argument defines the main title of the draw and group.names define the names of each group. If no output file and file extension is defined, the draw will be generated at your working directory as map.eps. The eps extension is only the default option but there are others that can be used. You can have access to a list of them on the help page.

We also can draw a map for a specific linkage group:

draw_map2(LG1_f2_final, col.group = "#58A4B0", col.mark = "#335C81", output = "map_LG1.pdf")

You can also change the default colors using the col.group and col.mark arguments.

With argument tag you can highlight some markers at the figure according to your specific purpose.

Speed up two-points ordering approaches

To use the parallelized approach you just need to include the arguments when using the functions:

# Without parallelization 
rcd_map <- rcd(input.seq = seq_all)

# With parallelization 
rcd_map_par <- rcd(input.seq = seq_all,
                   phase_cores = 4, 
                   size = batch_size, 
                   overlap = 30)

# Without parallelization 
record_map <- record(input.seq = seq_all)

# With parallelization 
record_map_par <- record(input.seq = seq_all,
                         phase_cores = 4, 
                         size = batch_size, 
                         overlap = 30)

# Without parallelization 
ug_map <- ug(input.seq = seq_all)

# With parallelization
ug_map_par <- ug(input.seq = seq_all,
                 phase_cores = 4, 
                 size = batch_size, 
                 overlap = 30)

Speed up mds_onemap

# Without parallelization ok
map_mds <- mds_onemap(input.seq = seq_all)

# With parallelization
map_mds_par <- mds_onemap(input.seq = seq_all, 
                          phase_cores = 4, 
                          size = batch_size, 
                          overlap = 30)

Speed up map estimation for an arbitrary order (map function)

Because we simulate this dataset we know the correct order. We can use map_overlapping_batches to estimate genetic distance in this case. This is equivalent to map, but with the parallelized process.

Similarly with map, using argument rm_unlinked = TRUE the function will return a vector with marker numbers without the problematic marker. To repeat the analysis removing automatically all problematic markers use map_avoid_unlinked:

# Without parallelization
batch_map <- map_avoid_unlinked(input.seq = seq_all)

# With parallelization
batch_map_par <- map_avoid_unlinked(input.seq = seq_all,
                                    size = batch_size,
                                    phase_cores = 4,
                                    overlap = 30)

As you can see in the above maps, heuristic ordering algorithms do not return an optimal order result, mainly if you don’t have many individuals in your population. Because of the erroneous order, generated map sizes are not close to the simulated size (100 cM) and their heatmaps don’t present the expected color pattern. Two of them get close to the color pattern, they are the ug and the MDS method. They present good global ordering but not local. If you have a reference genome, you can use its position information to rearrange local ordering. # Export estimated progeny haplotypes (new!)

Function progeny_haplotypes generates a data.frame with progeny phased haplotypes estimated by OneMap HMM. For progeny, the HMM results in probabilities for each possible genotype, then the generated data.frame contains all possible genotypes. If most_likely = TRUE the most likely genotype receives 1 and the rest 0 (if there are two most likely both receive 0.5), if most_likely = FALSE genotypes probabilities will be according to the HMM results. You can choose which individual to be evaluated in ind. The data.frame is composed of the information: individual (ind) and group (grp) ID, position in centimorgan (pos), progeny homologs (homologs), and from each parent the allele came (parents).

(progeny_haplot <- progeny_haplotypes(LG2_f2_final, most_likely = TRUE, ind = 2, group_names = "LG2_final"))

##    ind marker       grp      pos prob progeny.homologs parents allele
## 1 IND2    M23 LG2_final  0.00000    0               H1      P1  H1_P1
## 2 IND2    M60 LG2_final 20.56485    0               H1      P1  H1_P1
## 3 IND2    M29 LG2_final 23.27857    0               H1      P1  H1_P1
## 4 IND2   SNP1 LG2_final 31.78339    0               H1      P1  H1_P1
## 5 IND2    M44 LG2_final 34.44508    0               H1      P1  H1_P1
## 6 IND2   SNP2 LG2_final 36.34228    0               H1      P1  H1_P1
## 7 IND2    M36 LG2_final 39.16418    0               H1      P1  H1_P1
## 8 IND2    M40 LG2_final 46.11760    1               H1      P1  H1_P1
## 9 IND2    M26 LG2_final 49.76391    1               H1      P1  H1_P1
##  [ reached 'max' / getOption("max.print") -- omitted 135 rows ]

You can also have a view of progeny estimated haplotypes using plot. It shows which markers came from each parent’s homologs. position argument defines if haplotypes will be plotted by homologs (stack) or alleles (split). split option is a good way to better view the likelihoods of each allele.

plot(progeny_haplot, position = "stack")

plot(progeny_haplot, position = "split")