-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathGCLoess.R
64 lines (44 loc) · 1.71 KB
/
GCLoess.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
#!/urs/bin/env Rscript
####################################
#####Load required libraries########
####################################
library(dplyr)
library(ggplot2)
####################################
#####set file and path########
####################################
f <- "file.subset.GCcontent.txt"
setwd("/projects/dumont-lab/uma/k31/GCcontent")
path = "/projects/dumont-lab/uma/k31/GCcontent"
####################################
fname = file.path(path,f)
n = gsub("", "", f)
df = read.csv(fname, header = TRUE)
df <- df[,-1]
df[,2] <- as.numeric(df[,2])
loessMod30 <- loess(df[,2] ~ df[,4], data=df, span=0.3)
loessMod40 <- loess(df[,2] ~ df[,4], data=df, span=0.4)
smoothed30 <- predict(loessMod30)
smoothed40 <- predict(loessMod40)
smoothed30df <- cbind(smoothed30,df[,4])
smoothed30df <- as.data.frame(smoothed30df)
smoothed30df <- distinct(smoothed30df)
name <- paste(f,".subset1.smoothed30.txt", sep = "")
write.table(smoothed30df, file = name)
filename = paste(f,".subset1.Loess30.svg",sep="")
svg(filename = filename,
width = 9,
height = 6)
plot(smoothed30,x=df[,4],type="p",main="Loess Smoothing (30) and Prediction", xlab="GC content", ylab=f,ylim = c(0,50))
dev.off()
smoothed40df <- cbind(smoothed40,df[,4])
smoothed40df <- as.data.frame(smoothed40df)
smoothed40df <- distinct(smoothed40df)
name <- paste(f,".subset1.smoothed40.txt", sep = "")
write.table(smoothed40df, file = name)
filename = paste(f,".subset1.Loess40.svg",sep="")
svg(filename = filename,
width = 9,
height = 9)
plot(smoothed40,x=df[,4],type="p",main="Loess Smoothing (40) and Prediction", xlab="GC content", ylab=f,ylim = c(0,50))
dev.off()