calc.Tm.matrix <- function (sequences, conc=2e-07 ) { # this procedure calculates Tm for a sequences passed as a matrix of characters # K positions are in rows and probe sequences are in columns # N.B. this code assumes all sequences are equal length ## to handle next: missing bases are padded # "AA" "AT" "AC" "AG" "TA" "TT" "TC" "TG" "CA" "CT" "CC" "CG" "GA" "GT" "GC" "GG" dH.values <- c( -7.9, -7.2, -8.4, -7.8, -7.2, -7.9, -8.2, -8.5, -8.5, -7.8, -8.0, -10.6, -8.2, -8.4, -9.8, -8.0 ) # "AA" "AT" "AC" "AG" "TA" "TT" "TC" "TG" "CA" "CT" "CC" "CG" "GA" "GT" "GC" "GG" dS.values <- c( -22.2, -21.3, -22.4, -21.0, -21.3, -22.2, -22.2, -22.7, -22.5, -21.0, -19.9, -27.2, -22.2, -22.4, -24.4, -19.9 ) K = dim(sequences)[1] # number of bases in each sequence (or longest sequence) M = dim(sequences)[2] Bases = c("A","T","C","G") # order of bases for this function (not standard) dH <- rep( 0, M ) # Enthalpy vector dS <- rep( -10.8, M ) # Entropy seq.len <- apply( !is.na( sequences), 2, sum ) # lengths of individual sequences # computations for ends of probes dS <- dS + ifelse( sequences[1,] == "A" | sequences[1,]== "T", 4.1, -2.8 ) dH <- dH + ifelse( sequences[1,] == "A" | sequences[1,]== "T", 2.3, 0.1 ) dS <- dS + ifelse( sequences[K,] == "A" | sequences[K,]== "T", 4.1, -2.8 ) dH <- dH + ifelse( sequences[K,] == "A" | sequences[K,]== "T", 2.3, 0.1 ) # K by [seq.len,] == # test for palindrome - probably not necessary nn <- which( apply( sequences == apply( sequences, 2, rev ), 2, all ) ) dS[ nn ] <- dS[ nn ] - 1.4 cat("\nnumber of palindromes in probe sequences:",length(nn)) # now look at each base pair in interior of sequence and add for (kk in 1:(K-1)) { seq.index1 <- match( sequences[ kk , ], Bases) seq.index2 <- match( sequences[ kk+1, ], Bases) pair.index <- 4*(seq.index1 - 1) + seq.index2 dH <- dH + dH.values[ pair.index ] dS <- dS + dS.values[ pair.index ] } return( 1000*dH/(dS+1.987*log(conc/2)) - 273.15) }