`map.loc.lme` <-
function(input, s.chr, chrSet, prevLoc=NULL)
{
dfMerged <- input$dfMerged
fixed <- input$envModel$fixed
n.chr <- length(input$map)
map <- input$mapp[[s.chr]]
mrk <- grep(paste("C", s.chr, "M", sep=""), names(dfMerged))
chr <- sort(c(mrk, grep(paste("C", s.chr, "P", sep=""), names(dfMerged))))
type <- attr(input, "type")
results <- list()
formula <- list()
chrRE <- vector()
if (type=="f2") mrk <- mrk[seq(1, length(mrk), 2)]
wald <- rep(0, length(map))
f.pos <- vector()
f.mrk <- vector()
if (length(prevLoc)>0) {
f.pos <- prevLoc$pos
f.mrk <- prevLoc$mrk
if (type=="f2") {
f.pos <- paste(unique(substr(f.pos, 1, nchar(f.pos)-1)), "D", sep="")
f.mrk <- unique(substr(f.mrk, 1, nchar(f.mrk)-1))
}
}
# Set up random effects for markers on each chromosome
for (kk in 1:n.chr)
chrRE[kk] <- paste("pdIdent(~", paste(setdiff(names(dfMerged)[grep(paste("C", kk, "M", sep=""), names(dfMerged))], prevLoc$mrk), collapse="+"), "-1)", sep="")
int <- vector()
if (length(f.pos)>0) {
fmrk <- sapply(match(f.pos, names(dfMerged)), function(x) {
if (x==min(mrk)) return(c(min(mrk), min(mrk[mrk>x]))) else if (x==max(mrk)) return(c(max(mrk[mrk<x]), max(mrk))) else return(c(max(mrk[mrk<x]), min(mrk[mrk>x])))})
# because we want to exclude both additive and dominant effects
if (type=="f2") fmrk[2,] <- fmrk[2,]+1
int <- eval(parse(text=paste("c(", paste(apply(fmrk, 2, function(x) return(paste(x[1], ":", x[2], sep=""))), collapse=","), ")", sep="")))
}
# Loop over chrPos on the selected chromosome
for (jj in 1:length(map))
{
if (type=="f2") pos <- c(2*jj-1, 2*jj) else pos <- jj
if (all(!(chr[pos] %in% int)))
{
if (length(chrSet)>2)
formula$random <- paste("pdBlocked(list(", paste(chrRE[setdiff(chrSet,s.chr)], collapse=","), "))", sep="")
if (length(chrSet)==2)
formula$random <- chrRE[setdiff(chrSet, s.chr)]
formula$fixed <- paste(as.character(fixed)[2], "~", as.character(fixed)[3], sep="")
# If there are markers already mapped on other chromosomes, add in as fixed effects
if (length(f.pos) >0)
formula$fixed <- paste(formula$fixed, "+",paste(prevLoc$pos, collapse="+"), sep="")
effectnames <- names(dfMerged)[chr[pos]]
formula$fixed <- paste(formula$fixed, "+", paste(effectnames, collapse="+"), sep="")
formula$fixgrp <- paste(formula$fixed, "| grp1", sep="")
formula$fixgrp <- as.formula(formula$fixgrp)
formula$fixed <- as.formula(formula$fixed)
gd <- groupedData(formula$fixgrp, data=dfMerged)
# Fit model - different forms depending on relevant terms
if (length(chrSet)>1)
model <- lme(fixed=formula$fixed, random=eval(parse(text=formula$random)), data=gd, control=lmeControl(maxIter=input$maxit), na.action=na.omit)
if (length(chrSet)==1)
model <- lme(fixed=formula$fixed, random=~1|grp1, data=gd, control=lmeControl(maxIter=input$maxit), na.action=na.omit)
# Get output - Wald statistic for position fixed effect
############## need to put in a check for f2 here to get proper Wald####
wald[jj] <- anova(model, Terms=effectnames)[1,3]
} # end of check for distinct intervals
}
results$wald <- wald
return(results)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.