###################################################################################
#    Criar arquivos de LD (média do LD em janelas) à partir da Saída do Haploview #
###################################################################################

ld <- function (breed, nchromosomes)
{

# Abrir o pacote Hmisc para poder utilizar a função cut2
library('Hmisc')

for(i in 1:nchromosomes) {

# Importar os arquivos de saída do Haploview, excluindo a primeira linha (cabeçalho)
eval(parse(text=paste('chr',i,'<- read.table(\'',breed,'_CHR',i,'.haplo.LD\', skip=1, header=F)', sep='')))

# Criar linha de cabeçalho para cada arquivo no R
eval(parse(text=paste('names(chr',i,')<- c(\'L1\', \'L2\', \'D\', \'LOD\', \'r2\', \'CIlow\', \'CIhi\', \'Dist\', \'T-int\')', sep='')))

# Criar coluna com o número do cromossomo
eval(parse(text=paste('chr',i,'$chr <- ',i, sep='')))

# Remover os dados de LD de marcadores distantes mais de 10Mb
eval(parse(text=paste('chr',i,'<- subset(chr',i,', Dist <= 10000000)', sep='')))

# Criar classes de distância entre marcadores à partir de 0 até 10Mb, separados a cada 10Kb
eval(parse(text=paste('chr',i,'$distbin <- cut2(chr',i,'$Dist, cuts= seq(from=0, to=10000000, by=10000))', sep='')))

}

# Criar lista com os nomes dos arquivos que serão unidos com a função rbind
chr_list<- paste('chr', seq(1, nchromosomes), sep='', collapse=',')

# Unir os Dataframes de cada cromossomo para fazer um único arquivo
eval(parse(text=paste('data_distbin<- rbind(',chr_list,')', sep='')))

# Gerar um novo arquivo com médias e desvios-padrão do desequilíbrio de ligação par cada classe de distância entre os marcadores
data_distbin2<- data.frame(distbin=tapply(as.character(data_distbin$distbin), as.factor(data_distbin$distbin), head, n=1, ordered=T), n_obs=tapply(data_distbin$D, data_distbin$distbin, length), mean_D=tapply(data_distbin$D, data_distbin$distbin, mean), sd_D=tapply(data_distbin$D, data_distbin$distbin, sd), mean_r2=tapply(data_distbin$r2, data_distbin$distbin, mean), sd_r2=tapply(data_distbin$r2, data_distbin$distbin, sd))

# Remover as linhas contendo NAs
data_distbin2<- data_distbin2[complete.cases(data_distbin2),]

# Gerar o arquivo de saída contendo o dataframe data_distbin
write.table(data_distbin, file= paste(breed, '_LD_allchromosomes.txt', sep=''), quote=F, sep='\t', row.names=F)

# Gerar o arquivo de saída contendo o dataframe data_distbin2
write.table(data_distbin2, file= paste(breed, '_LD_MeanByDistbin.txt', sep=''), quote=F, sep='\t', row.names=F)

}