################################################### ### chunk number 1: ################################################### library(GeneR) s <- "gtcatgcatgctaggtgacagttaaaatgcgtctaggtgacagtctaacaa" s2 <- insertSeq(s,"aaaaaaaaaaaaaa",20) s2 ################################################### ### chunk number 2: ################################################### library(GeneR) strComp(s) ################################################### ### chunk number 3: ################################################### strCompoSeq(s2,wsize=1) strCompoSeq(s2,wsize=2) ################################################### ### chunk number 4: ################################################### seqNcbi("BY608190",file='toto.fa') strReadFasta('toto.fa',from=10,to=35) ################################################### ### chunk number 5: ################################################### readFasta('toto.fa') readFasta('toto.fa',seqno=1) ################################################### ### chunk number 6: ################################################### getSeq(from=10,to=35) ################################################### ### chunk number 7: ################################################### getSeq(seqno=0,from=c(1,10,360),to=c(10,20,0)) assemble(seqno=0,destSeqno=1,from=c(1,10,360),to=c(10,20,0)) getSeq(seqno=1) ################################################### ### chunk number 8: ################################################### sizeSeq(seqno=0) size0 <- sizeSeq(seqno=0) getSeq(seqno=0,from=size0-20,to=0) ################################################### ### chunk number 9: ################################################### getSeq(0,from=1,to=35) getSeq(1) appendSeq(0,1) getSeq(seqno=0,from=size0-20,to=size0+10) concat(seqno1=0,seqno2=1,destSeqno=3, from1=2,to1=10,from2=8,to2=0, strand1=1) getSeq(3) ################################################### ### chunk number 10: ################################################### exactWord("AAAG",seqno=0) ################################################### ### chunk number 11: ################################################### getOrfs(seqno=0) maxOrf(seqno=0) ################################################### ### chunk number 12: ################################################### x=c(5,15,30);y=c(10,20,35) lowerSeq(from=x,to=y) getSeq(from=1,to=35) ################################################### ### chunk number 13: ################################################### posMaskedSeq(seqno=0,type="lower") posMaskedSeq(seqno=0,type="upper") ################################################### ### chunk number 14: ################################################### s <- "ATGCtgTGTTagtacATNNNNNNNNNNNNNNNTGGGTTTaAAAattt" placeString(s,upper=FALSE,seqno=0) posMaskedSeq(seqno=0,type="N") ################################################### ### chunk number 15: ################################################### dnaToRna() getSeq() ################################################### ### chunk number 16: ################################################### s <- "xxxxxxxxxxATGTGTCGTTAATTGxxxxxxxxxxxxxxx" placeString(s) setStrand(0) writeFasta(file="toto.fa") readFasta(file="toto.fa",from=11,to=25) getSeq() ################################################### ### chunk number 17: ################################################### getOrfs() AtoT(getOrfs()) ################################################### ### chunk number 18: ################################################### getSeq() exactWord(word="AA") AtoT(exactWord(word="AA")[[1]]) setStrand(1) getSeq() exactWord(word="AA") AtoT(exactWord(word="AA")[[1]]) exactWord(word="TT") AtoT(exactWord(word="TT")[[1]]) ################################################### ### chunk number 19: eval=FALSE ################################################### ## seqNcbi("BY608190",file="bank.fa") ## seqNcbi("BY608190",file="BY608190.gbk",type="genbank") ## ## #seqSrs("embl|BY608190|BY608190",file="BY608190.embl",submotif=TRUE) ################################################### ### chunk number 20: ################################################### setStrand(0) seqNcbi("BY608190",file="bank.fa") seqNcbi("BY608190",file="BY608190.gbk",type="genbank") #seqSrs("embl|BY608190|BY608190",file="BY608190.embl",submotif=TRUE) ################################################### ### chunk number 21: ################################################### readFasta(file="bank.fa", name="gi|26943372|gb|BY608190.1|BY608190", from=1,to=30) getSeq() fastaDescription(file="bank.fa", name="gi|26943372|gb|BY608190.1|BY608190") ################################################### ### chunk number 22: ################################################### s="cgtagctagctagctagctagctagctagcta" placeString(s,seqno=0) writeFasta(seqno=0,file="bank.fa",name="MySequence", comment="A sequence Generated by R",append=TRUE) ## Show bank file cat(paste(readLines("bank.fa"),collapse='\n')) ################################################### ### chunk number 23: ################################################### # readEmbl(file='BY608190.embl',name='BY608190',from=10,to=30) # getSeq() s<-"gtcatgcatcctaggtgtcagggaaaatgcgtctacgtgacagtctaacaa" placeString(s) # Add lines with "CC bla bla bla" and a line "XX" writeEmblComment(file="toto.embl",code="CC",text="This is a comment for \ this dummy sequence... I try to be long enough to show that this comment \ will be written on several lines",append=FALSE) # Add a line with "FT CDS bla bla bla" writeEmblLine(file="toto.embl",code="FT",header="CDS",text="<1..12", nextfield = FALSE) # Add lines with "FT bla bla bla" writeEmblLine(file="toto.embl",code="FT",header="",text="/codon_start=2", nextfield = FALSE) writeEmblLine(file="toto.embl",code="FT",header="",text="/gene=\"toto\"", nextfield = FALSE) writeEmblLine(file="toto.embl",code="FT",header="",text="/note=\"Here is \ what I think about this gene\"",nextfield = TRUE) ## Translation prot <- translate(seqno=0,from=getOrfs()[1,1],to=getOrfs()[1,2]) writeEmblLine (file="toto.embl",code='FT',header='', text=paste('/translation="',prot ,'\"', sep=''),nextfield =FALSE) # Add sequence writeEmblSeq(file="toto.embl") ## Show file cat(paste(readLines("toto.embl"),collapse='\n')) ################################################### ### chunk number 24: ################################################### readGbk(file='BY608190.gbk',name='BY608190',from=10,to=30) getSeq() ################################################### ### chunk number 25: ################################################### a = matrix(c(1,5,15,45,17,38, 100,120,130,140, 135,145,142,160), ncol=2,byrow=TRUE) b = matrix(c(15,18, 28,45, 1,10, 15,20, 25,40, 17,23, 35,38,100,105, 110,120),ncol=2,byrow=TRUE) a <- as.segSet(a) b <- as.segSet(b) A = list( matrix( c( 1, 15, 17, 5, 45, 38),ncol=2), matrix( c( 100 , 120),ncol=2), matrix( c( 130, 135, 140, 145),ncol=2), matrix( c( 142 , 160),ncol=2)) B = list( matrix( c(15, 28, 18, 45),ncol=2), matrix( c(1, 15, 25, 10, 20, 40),ncol=2), matrix( c(17, 35, 23, 38),ncol=2), matrix( c(100, 110, 105, 120),ncol=2)) A = as.globalSeg(A) B = as.globalSeg(B) c = or(a,b) d = and(a,b) e = not(a,b) f = Xor(a,b) par(mfrow=c(4,1)) ##par(mfrow=c(2,1)) plot(a,xlim=c(1,160),main="a (segSet)") plot(b,xlim=c(1,160),main="b (segSet)") plot(A,xlim=c(1,160),main="A (globalSet)") plot(B,xlim=c(1,160),main="B (globalSet)") ################################################### ### chunk number 26: ################################################### g = or(a) h = and(a) #i = not(a) #j = Xor(a) C = or(A,B) D = and(A,B) E = not(A,B) F0 = Xor(A,B) G = or(A) H = and(A) #I = not(A) #J = Xor(A) par(mfrow=c(4,1)) ##par(mfrow=c(2,1)) plot(c,xlim=c(1,160),main="a or b") plot(d,xlim=c(1,160),main="a and b") plot(e,xlim=c(1,160),main="a not b") plot(f,xlim=c(1,160),main="a Xor b") ################################################### ### chunk number 27: ################################################### par(mfrow=c(4,1)) ##par(mfrow=c(2,1)) plot(C,xlim=c(1,160),main="A or B") plot(D,xlim=c(1,160),main="A and B") plot(E,xlim=c(1,160),main="A not B") plot(F0,xlim=c(1,160),main="A Xor B") ################################################### ### chunk number 28: ################################################### par(mfrow=c(4,1)) ##par(mfrow=c(2,1)) plot(G,xlim=c(1,160),main="or A") plot(H,xlim=c(1,160),main="and A") #plot(I,xlim=c(1,160),main="not A") #plot(J,xlim=c(1,160),main="Xor A") plot(g,xlim=c(1,160),main="or a") ##plot.segSet(h,xlim=c(1,160),main="and a") #plot(i,xlim=c(1,160),main="not a") #plot(j,xlim=c(1,160),main="Xor a") ################################################### ### chunk number 29: ################################################### ## Set seed of your choice (not requiered) set.seed(3) #### ---- RANDOMSEQ ---- ## Create a sequence of size 30, GC rich randomSeq(prob = c(0.20, 0.30, 0.20, 0.30), letters = c("T", "C","A", "G"), n = 30) ## [1] "CTGGAACCGAGGGGTTCATCCCCCCAGTGA" ## use with bi-nucleotides randomSeq(prob=rep(0.0625,16),letters = c("TT","TC","TA","TG", "CT","CC","CA","CG","AT","AC","AA","AG","GT","GC","GA","GG"),n=10) ## [1] "CGCATGATCCCAGGCTAACT" #### ---- SHUFFLESEQ ---- ## Create a sequence with 7 T, 3 C and A, and 4 G. shuffleSeq(count=c(7,3,3,4,0),letters=c("T","C","A","G","N")) ## [1] "TATCTTTTGTCGGACGA" ## Same with bi-nucleotides shuffleSeq(count=c(rep(4,4),rep(2,4),rep(1,4),rep(0,4)), letters = c("TT","TC","TA","TG","CT","CC","CA", "CG","AT","AC","AA","AG","GT","GC","GA","GG")) ## [1] "TCTTTCCATTCCTTCTAGTGTACCCGTATACGTGTCTGTGTACTTCAACAACTTAT" ## From a real sequence: seqNcbi("BY608190",file="BY608190.fa") readFasta("BY608190.fa") ## create a random sequence from a real tri-nucleotides distribution ## Size of sequence will be 10*3. randomSeq(compoSeq(wsize=3,p=TRUE),n=10) ## re assemble real tri-nucleotides of a real sequence shuffleSeq(compoSeq(wsize=3,from=1,to=30,p=FALSE)) ################################################### ### chunk number 30: ################################################### ## We create a sequence with a biais every 100 nucleotides s <- "" for(i in 1:20) s <- paste(s, randomSeq(n=100),randomSeq(prob=c(0.3,1,1,1,0)/3.3,n=100),sep="") placeString(s,seqno=0) dens <- densityProfile(ori=200*(1:20)-100,from=200*(0:19)+1,to=200*(1:20), seqno=0, fun=seqSkew,nbinL=10,nbinR=10,sizeBin=10) plot(dens$skta,main="TA skew") ## Example with flagged 'N' ## We create a sequence with a biais every 100 nucleotides s <- "" for(i in 1:20) s <- paste(s, randomSeq(n=100),randomSeq(prob=c(0.3,1,1,1,0.2)/3.5,n=100),sep="") placeString(s,seqno=1) dens2 <- densityProfile(ori=200*(1:20)-100,from=200*(0:19)+1,to=200*(1:20), seqno=1, fun=compoSeq,nbinL=10,nbinR=10,sizeBin=10) plot(dens2$T,main="#T") ## The same but more permissive (allow 4 N in each bin) dens3 <- densityProfile(ori=200*(1:20)-100,from=200*(0:19)+1,to=200*(1:20), seqno=1, fun=compoSeq,nbinL=10,nbinR=10,sizeBin=10,threshold=4) plot(dens3$T,main="#T")