view from_imgt.r @ 0:5560672b1ca4 draft default tip

Uploaded
author davidvanzessen
date Fri, 24 Jul 2015 04:44:39 -0400
parents
children
line wrap: on
line source

library(data.table)

args <- commandArgs(trailingOnly = TRUE)

infile="D:/wd/prisca/Mouse data Groningen July 2015/JIVFXVQ01_MAAIKE_1_PB_IGH_MID8_10nt_trimmed/1_Summary.txt"
patient="JIVFXVQ01"
sample="sample1"
cell.count=10000
receptor="IgH"
output="D:/wd/prisca/mousetest.txt"

infile=args[1]
patient=args[2]
sample=args[3]
cell.count=args[4]
receptor=args[5]
output=args[6]

dat = read.table(infile, header=T, sep="\t", fill=T, stringsAsFactors=F)
dat = dat[,c("V.GENE.and.allele", "J.GENE.and.allele", "AA.JUNCTION", "Sequence")]

dat = dat[dat$V.GENE.and.allele != "",]
dat = dat[dat$J.GENE.and.allele != "",]
dat = dat[dat$Sequence != "",]

dat$V.GENE.and.allele = as.factor(as.character(lapply(strsplit(as.character(dat$V.GENE.and.allele), ", "), "[[", 1)))
dat$J.GENE.and.allele = as.factor(as.character(lapply(strsplit(as.character(dat$J.GENE.and.allele), ", "), "[[", 1)))

dat$V.GENE.and.allele = gsub("Homsap ", "", dat$V.GENE.and.allele)
dat$V.GENE.and.allele = gsub("\\*.*", "", dat$V.GENE.and.allele)

dat$J.GENE.and.allele = gsub("Homsap ", "", dat$J.GENE.and.allele)
dat$J.GENE.and.allele = gsub("\\*.*", "", dat$J.GENE.and.allele)

dat = data.frame(data.table(dat)[, list(Clone_Molecule_Count_From_Spikes=.N), by=c("V.GENE.and.allele", "J.GENE.and.allele", "AA.JUNCTION", "Sequence")])

dat = dat[order(-dat$Clone_Molecule_Count_From_Spikes),]
dat$perc = 100 / nrow(dat) * dat$Clone_Molecule_Count_From_Spikes

dat$Log10_Frequency = log10(dat$perc / 100)

dat$Patient = patient
dat$Sample = sample
dat$Receptor = receptor
dat$Cell_Count = cell.count
dat$Total_Read_Count = dat$Clone_Molecule_Count_From_Spikes
dat$Related_to_leukemia_clone = F

dat = dat[,c("Patient", "Receptor", "Sample", "Cell_Count", "Clone_Molecule_Count_From_Spikes", "Log10_Frequency", "Total_Read_Count", "V.GENE.and.allele", "J.GENE.and.allele", "Sequence" ,"AA.JUNCTION", "Related_to_leukemia_clone")]

names(dat) = c("Patient", "Receptor", "Sample", "Cell_Count", "Clone_Molecule_Count_From_Spikes", "Log10_Frequency", "Total_Read_Count", "V_Segment_Major_Gene", "J_Segment_Major_Gene", "Clone_Sequence" ,"CDR3_Sense_Sequence", "Related_to_leukemia_clone")

write.table(dat, output, quote=F, sep="\t", na="", dec=".", row.names=F, col.names=F)

output