Reputation: 13
I am attempting to subset a region of genetic data by its chromosome and position information. Unfortunately, my results are not within my parameters. Any help would be appreciated.
Here is my code:
subset.by.region<- function(df,region.info, expansion=0){
MBstart = as.numeric(region.info[[3]]) - expansion
MBend = as.numeric(region.info[[4]]) + expansion
chrom =as.numeric(region.info[[2]])
print(chrom)
print(MBstart)
print(MBend)
BPstart <- MBstart * 1e6
BPend <- MBend * 1e6
sub_results <- as.numeric(df$CHR) == chrom & as.numeric(df$BP) >= BPstart & as.numeric(df$BP) <= BPend;
print(head(sub_results))
region_results <- subset(results, sub_results)
return(region_results)
}
Here are the prints to the console containing the region information being used:
[1] 1
[1] 113.308
[1] 115.158
Here is the print from the subset (region_results):
GENE CHR SNP EMP1 NP BP SNP_IM SNP_LZ
3238 AP3S1 5 rs26538 1.00000 6 115178395 rs26538 rs26538
3239 AP4B1 1 rs1217401 1.00000 46 114438951 imm_1_114240474 rs1217401
3240 AP4B1 1 rs1217402 1.00000 41 114440258 imm_1_114241781 rs1217402
3241 AP4B1 1 rs3789613 1.00000 297 114443035 imm_1_114244558 rs3789613
3242 AP4B1 1 rs7523862 1.00000 297 114443419 imm_1_114244942 rs7523862
3243 AP4B1 1 rs17464525 1.00000 148 114443899 imm_1_114245422 rs17464525
As you can see, there is a row in the subset containing a marker in chromosome 5. What am I doing wrong? Thank you in advance. edit: Here is the call to the function with the stuff before it:
write.genelist <- function(table_loc, region.info, out_folder,yank_loc){
region.ID = as.character(region.info[[1]])
out_name = paste0(region.ID,"_genes.list")
region_folder = file.path(out_folder, region.ID)
out_loc <- file.path(region_folder,out_name, fsep = .Platform$file.sep)
results <- read.table(table_loc, T,strip.white = TRUE)
gene_region_results <- subset.by.region(results,region.info)
...
}
Upvotes: 1
Views: 706
Reputation: 174788
I would use [
to subset in a function not subset()
. See ?subset
for why.
subset.by.region<- function(results, df, region.info, expansion=0){
MBstart = as.numeric(region.info[[3]]) - expansion
MBend = as.numeric(region.info[[4]]) + expansion
chrom =as.numeric(region.info[[2]])
print(chrom)
print(MBstart)
print(MBend)
BPstart <- MBstart * 1e6
BPend <- MBend * 1e6
sub_results <- as.numeric(df$CHR) == chrom &
as.numeric(df$BP) >= BPstart & as.numeric(df$BP) <= BPend
print(head(sub_results))
results[sub_results, ]
}
I would also pass in results
not rely on it being found in the global environment.
Upvotes: 1