Question: Loop in R script
0
gravatar for Seq225
11 days ago by
Seq22590
Seq22590 wrote:

Hello,

I have this R script:

#!/usr/bin/env Rscript

library(ggplot2)

Pai <- read.csv("abc.txt", sep = "\t")

colnames(Pai) <- c("nt","a","c","g","t")

Pai[2:5] <- scale(Pai[2:5])

write.table(Pai, file = "PaiZscore.txt", row.names = FALSE)

pdf("abc.pdf")

ggplot() + 
    geom_line(aes(Pai$nt,Pai$t), color = 'red') + 

theme_classic()

dev.off()

It runs one file abc.txt and create one pdf abc.pdf. How can I modify this script so that it can run as a loop, use 200 files as input (abs1.txt, abc2.txt.........abc200.txt) and creates 200 pdfs (abs1.pdf, abc2.pdf.........abc200.pdf)?

The abc.txt files is generated using multiple custom pipelines, which incorporated bowtie, samtools, bedtools, and piPipes. I can use the aforementioned tools, but I am no expert in R. My goal is to make the heatmap of the final result in R.

Thanks a ton!!

ADD COMMENTlink modified 10 days ago by bioinformatics2020330 • written 11 days ago by Seq22590
1

Check For and While loops in R https://cran.r-project.org/doc/contrib/Paradis-rdebuts_en.pdf

ADD REPLYlink written 11 days ago by JC11k
1
  1. write a function for ggplot
  2. Loop the function over the text files
ADD REPLYlink written 11 days ago by cpad011213k

Or,

  1. Write a function that accepts a file name and outputs a ggplot
  2. lapply that function over a list of file names. Avoid loops.
ADD REPLYlink written 11 days ago by RamRS28k
3
gravatar for bioinformatics2020
10 days ago by
bioinformatics2020330 wrote:

Put this in the directory where you have the 200 files. Make sure you don't have any other text files because it uses the list.files() function to look for .txt files.

if(!require("tidyverse")) install.packages("tidyverse")
if(!require("magrittr")) install.packages("magrittr")
library(tidyverse) 
library(magrittr)
files <- list.files(pattern = "*\\.txt$")
Pai <- map(files, read_tsv)
Pai <- map(Pai, ~ bind_cols(.x[1], as_tibble(scale(.x[2:5]))))
iwalk(Pai, ~ write_tsv(.x, path = paste0("zscore_",.y,files[[.y]])))

plot.me <- function(x, x.cord, y.cord, filename, device = "pdf", height = 10, width = 10) { 
  ggplot(data = x, mapping = aes_string(x = x.cord, y = y.cord)) + 
    geom_point(color = "red") +
    theme_classic()
  ggsave(filename = paste(filename, device, sep = "."), device = device, height = height, width = width)
}
iwalk(Pai, ~plot.me(.x, x.cord = "nt", y.cord = "t", device = "pdf", filename = gsub(".txt","",files)[.y]))
ADD COMMENTlink modified 8 days ago • written 10 days ago by bioinformatics2020330

Interesting. Why use map instead of lapply?

ADD REPLYlink written 10 days ago by RamRS28k
1

The main reason here was because I wanted to use imap, which allows you to iterate over the list of objects and their indices, and so I just wanted to be consistent throughout the calls.

ADD REPLYlink written 10 days ago by bioinformatics2020330
1

Thanks for the explanation - I've never heard of imap before - learned something new today!

ADD REPLYlink written 10 days ago by RamRS28k
1

If you want to use imap to save a file, or any other purpose where you don't want to return a list, there's a handy alternative function called iwalk. There is also walk, pwalk, and walk2 that are alternatives to map, pmap, and map2.

ADD REPLYlink written 10 days ago by rpolicastro720

Thanks a lot!!

I am no expert in R. Even understanding a minor issue takes time for people like me.

Would you please help me fix this or direct me to somewhere where I can get some helpful insights? I have tried, but it seems to me that x needs to be defined as a function??

Error in .x[1] : object of type 'closure' is not subsettable
Calls: map -> .f -> bind_cols -> list2
Execution halted

Thanks!

ADD REPLYlink modified 8 days ago • written 9 days ago by Seq22590
1

Can you send me what one of these files look like inside?

ADD REPLYlink written 8 days ago by bioinformatics2020330

Thanks for your time, and appreciate your help. I don't quite see the option of sending you the file. I am pasting it below.

Thanks!!

"nt" "a" "c" "g" "t"
2 1.25468332957509 0.158563316582844 -1.64006710598786 0.926245535065208
3 0.219803815359422 -0.514600017160217 -1.56107070301328 2.21008499426014
4 -0.610531627769753 -0.333275269517582 -0.450825748677823 1.3783121894941
5 -0.562435399693933 -1.41509014868781 1.52408432568676 -0.593822373913796
6 -0.936667286688371 0.143229299868879 1.08652268935237 -0.677621308205621
7 -0.390315752029343 0.86009458124674 0.534481632158826 -0.910368496624458
8 -0.341949320425008 0.678386483186255 -0.0672356500731077 -0.0358259093712009
9 0.472984521578764 0.0979939505626822 -0.443542392375202 0.0917374535121999
10 1.15281659932462 -0.368543507959701 0.189922853227132 -1.07473376536896
11 2.53625866532571 -0.486615436657231 -1.05366354726401 -0.615455926683496
12 1.06067719610071 -0.278456159765157 -0.794077258529559 0.261821937356734
13 0.215210355374653 -0.700524969817042 -0.448211210517908 0.853139046395189
14 0.836408267432568 -0.226703853355526 -0.50610455548746 0.0512056247597743
15 0.447585389898275 -0.11974908677562 -0.850289828967738 0.797936187603542
16 0.527025227281932 0.235233400152669 -0.563064136828472 0.112127698651457
17 1.31385790232006 0.581398827470427 -0.987926587814707 -0.270811051524834
18 -0.0979555341752051 -0.136233154743132 -0.784552869518439 1.22314739721488
19 1.66323106469104 -0.89411693083085 0.0464967598832076 -1.01256838384683
20 0.106588536911285 -1.38633886734912 1.10146290740903 -0.665436893427285
21 0.546479881335072 1.30976462138376 -0.92947298466803 -0.114900274667344
22 -1.05906948510605 0.983533415794158 0.156867620776775 0.127793374795032
23 -1.34251298651916 -1.34608707347497 2.4569141905994 -1.16275994560429
24 0.269251061077821 -0.00321055974948641 1.26076298243815 -1.92441020001336
25 0.984209597530732 -1.67998528742155 -0.24054217953035 0.504266925293022
26 -1.21740875281633 -0.451347198215112 0.209718642152205 1.13387790934911
27 -1.15877458712839 3.22996686439002 -1.76537818493809 1.32186602307201
28 0.690228158505499 0.590599237498806 0.673612412811462 -1.91520972354809
29 0.0282295136416919 -0.948552690165425 1.68375790616731 -1.65262315199863
30 -1.45194541556808 -0.681357448924586 1.6882399715843 -0.469740272395635
31 -1.50463510362867 1.45237097682364 -0.253988375781343 0.780778542303435
32 0.900176300162306 1.43358680634903 -1.10464704138235 -0.287471373772764
33 -1.09419594381311 1.61874505817016 -0.997637729551536 1.285312778737
34 -1.10878693435296 -0.463614411586284 -0.493405370139301 1.97808379041957
35 0.855592717957192 -1.40895654200222 0.93413246517445 -1.11725488633009
36 0.228180124743413 -0.163067683992571 -0.186944147252086 0.144702358569051
37 -1.21092386813195 0.728222037506642 0.167512526142144 0.418976021844667
38 -1.46896823786457 0.0255407215891979 1.0826008821125 -0.106197121254247
39 0.385168374811116 0.768473831380799 -0.0661151337188583 -0.764901503862685
40 -0.965849267768082 0.273568441937581 0.944030359636986 -0.545582037852627
41 -0.171721097460029 -1.16284557374309 0.447641614704497 0.316278811570116
ADD REPLYlink modified 8 days ago by RamRS28k • written 8 days ago by Seq22590
1

I made some minor adjustments! Please try now. :)

ADD REPLYlink written 8 days ago by bioinformatics2020330

Thank you so much! It works now!!! Much appreciation for your help!!

ADD REPLYlink written 7 days ago by Seq22590
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1615 users visited in the last hour