-
Notifications
You must be signed in to change notification settings - Fork 0
/
gam.R
50 lines (39 loc) · 1.65 KB
/
gam.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
#!/home/ljw/miniconda3/envs/py310/bin/R
# -*- coding: utf-8 -*-
# This file is part of PCR efficiency calculator.
# PCR efficiency calculator is free software: you can
# redistribute it and/or modify it under the terms of the GNU General Public
# License as published by the Free Software Foundation, version 2.
#
# This program is distributed in the hope that it will be useful, but WITHOUT
# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
# details.
#
# You should have received a copy of the GNU General Public License along with
# this program; if not, write to the Free Software Foundation, Inc., 51
# Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
#
# Copyright Izaskun Mallona
# izaskun.mallona@gmail.com
# loading the library
suppressMessages(library(mgcv))
# loading the model
load('gam.RData')
dat<-read.table('./primerDataForR.dat',sep='\t',header=F)
# this is based on
#
# miniGaussian<-(gam(efficiency ~ s(lengthSequence,gcSequence)+s(primersLength,gcPrimers)+s(gcImbalance,primerDimers), data=dataMini))
#
# a = query.getLengthSequence()
# b = query.gcContent()
# c = query.getPrimersLength()
# d = query.getGcPrimers()
# e = query.getGcImbalance()
# f = query.getPrimerDimers()
colnames(dat)<-c('lengthSequence','gcSequence','primersLength','gcPrimers', 'gcImbalance','primerDimers')
efficiencies <- data.frame()
efficiencies <- predict(miniGaussian,dat)
# efficiencies[as.numeric(efficiencies)<1.5]<-'1.5'
# efficiencies[as.numeric(efficiencies)>2]<-'2'
write.table(efficiencies, './gamResult.data', row.names = FALSE,quote=FALSE)