| Author | Message | ||
     Howard Mark (hlmark) Senior Member Username: hlmark Post Number: 412 Registered: 9-2001 |
I don't know if this is coincidence or not, but I recently received an e-mail message offering a course in R programming. The message is promoting a course entitiled "The Fundamental of using R" and includes links to more infornation. The full e-mail is rather long to copy and post it all here, but if anyone is interested, they can send me a message asking for it, and I'll forward the e-mail. I have no idea whether the message is legitimate, or is spam or some type of on-line malware (virus, phishing, etc, so you select their links at your own risk (although my copy of MacAfee didn't raise any red flags about their web site.)) Howard \o/ /_\ | ||
| Magdalena Sut (vickym) Junior Member Username: vickym Post Number: 6 Registered: 11-2010 |
I am sorry, I think I am not long enough in R environment to do it...could someone please tell me, if I hav a table of 100 colums and 106 rows, how would the commend look like (lets say i previously uploded the table to R). what do I put in T? how i choose filter length? | ||
| Magdalena Sut (vickym) New member Username: vickym Post Number: 5 Registered: 11-2010 |
Thank you very much for ur help, now I will try to proceed.. Magda | ||
| Christian Mora (cmora) Senior Member Username: cmora Post Number: 26 Registered: 2-2007 |
Hi Magdalena Not sure to understand your problem. If you want to get the transformed spectra then the R code posted by Peter is fine. You can also type on the R environment ?sav_gol to display the help page for the function and an example of it. Another option is to post your question directly on the R-help list. BTW: you don't need the load the Rcmdr library to read the data into R. Now, if you want to understand what's behind the code of the sav_gol function then I guess you'll need to be familiar with the R language. Despite the obvious differences between Matlab and R, you can easily translate one code into the other. For those interested there is a nice publication on: http://www.math.umaine.edu/~hiebeler/comp/matlabR.html Regards CM | ||
| David W. Hopkins (dhopkins) Senior Member Username: dhopkins Post Number: 187 Registered: 10-2002 |
Hi Venky, Thanks for your comments. I will be interested in seeing what your R code looks like. Unfortunately, a lot of Peter's code was not useful for me. Best regards, Dave | ||
| venkatarman (venkynir) Senior Member Username: venkynir Post Number: 131 Registered: 3-2004 |
David ; R is GLP package and it is structure oriented . Matlab is function oriented . R-is flexible that means it works under Labview ,Vb and other environment . There is built in function in R. I will mail R-code soon. you can write smoothening using S.Golay publication. | ||
| David W. Hopkins (dhopkins) Senior Member Username: dhopkins Post Number: 186 Registered: 10-2002 |
Hi Magda, Oh my goodness! I see from Peter's post, that R doesn't look very much like Matlab. Best wishes, Dave | ||
| David W. Hopkins (dhopkins) Senior Member Username: dhopkins Post Number: 185 Registered: 10-2002 |
Hi Magda, I haven't used R, but I think it is a language similar in structure to Matlab, which I have used extensively and have implemented the Savitzky-Golay method. It is a 2-step procedure in Matlab. First, you have to obtain the convolution coefficients, by specifying the order of the 'derivative', 0 for smoothing, 1 for first derivative, 2 for second derivative, etc. You must also specify the order of the polynomial used to fit the data, and this is usually 1 for linear, 2 for quadratic, 3 for cubic, etc. The third parameter is the number of points used in the fitting polynomial, which is normally an odd integer such as 5, 7, 9.... The algorithm for obtaining the coefficients is not something that you should necessarily be able to follow in the R code. However, you should have some help that explains the 3 input integers that are required, and the output vector that is generated. The call should be something like: p = savgol(a,b,c) where a, b, c are the 3 parameters I've listed above, and p is the vector of coefficients. The process of convolution to obtain your treated spectrum is probably something that you will have to code for yourself. You need to take care that for a convolution of np points, there are (np-1)/2 spectral points at the beginning and end of the spectra that are undefined. I usually prefer to set these to zero, to emphasize that the values are not defined. However, in the case of a smooth, you may prefer to set these regions to the original data points, to minimize the discontinuity when the spectra are displayed at the original length of the spectra. Then you have to generate the convoluted spectrum point by point over the defined wavelength region by performing the multiplication and summing of terms over the proper interval. This is probably most easily done using vector multiplication. If your spectrum, s, is a row vector, then your convolution vector should be a column vector, the result r can be calculated fairly easily. I assume you will be able to interpret the Matlab code and make something equivalent in R. In Matlab, the notation a:b means a vector with the values from a to b at successive integral values. Code for a single spectrum that is npts long, using a convolution function of np points, would resemble this: r = zeros(npts); nz = (np-1) ./ 2; for i = 1 : npts - np +1 r(i+nz) = s(i : i+np-1) * p ; end The first line sets up the result vector of the proper length and zeroes out the ends (and in fact the rest of the vector, which will be replaced by the proper values). The FOR loop does the point-by-point calculation of the convoluted spectrum. You can check that you are obtaining the correct p coefficients by comparing the tables for smoothing, 1st and 2nd derivatives in the original S-G publication for quadratic fits of 5 points, but be warned that there are numerous tables that are erroneous. I know this is somewhat long, but I wanted to indicate to you that there could be a reason you could not understand the code of the part you have obtained in RTisean, and that you may need to do your own coding in order to obtain the desired results. I hope this provides just the hints you need. If you want to treat a number of spectra, I hope it will be obvious to you how to generalize the calculations sketched above. Finally, I wish you great success in your studies. I have had a great deal of fun getting into the use and understanding of Savitzky-Golay convolution functions, and I hope you will, too. Please feel free to ask further questions, if I have not been clear enough or missed your need. Perhaps a user of R can provide better code for you. Best wishes, Dave | ||
| Peter Tillmann (tillmann) Intermediate Member Username: tillmann Post Number: 19 Registered: 11-2001 |
Magdalena, here is what I use. I found it somewhere in the nabble R archive. # ---------------------------------------------------------------------- # Savitzky-Golay Algorithm # ---------------------------------------------------------------------- # T2 <- sav.gol(T, fl, forder=4, dorder=0); # # Polynomial filtering method of Savitzky and Golay # See Numerical Recipes, 1992, Chapter 14.8, for details. # # T = vector of signals to be filtered # (the derivative is calculated for each ROW) # fl = filter length (for instance fl = 51..151) # forder = filter order (2 = quadratic filter, 4= quartic) # dorder = derivative order (0 = smoothing, 1 = first derivative, etc.) # sav.gol <- function(T, fl, forder=4, dorder=0) { m <- length(T) dorder <- dorder + 1 # -- calculate filter coefficients -- fc <- (fl-1)/2 # index: window left and right X <- outer(-fc:fc, 0:forder, FUN="^") # polynomial terms and coefficients Y <- pinv(X); # pseudoinverse # -- filter via convolution and take care of the end points -- T2 <- convolve(T, rev(Y[dorder,]), type="o") # convolve(...) T2 <- T2[(fc+1) : (length(T2)-fc)] } #----------------------------------------------------------------------- # *** PseudoInvers of a Matrix *** # using singular value decomposition # pinv <- function (A) { s <- svd(A) # D <- diag(s$d); Dinv <- diag(1/s$d) # U <- s$u; V <- s$v # A = U D V' # X = V Dinv U' s$v %*% diag(1/s$d) %*% t(s$u) } #-------------------------- Peter | ||
| Magdalena Sut (vickym) New member Username: vickym Post Number: 4 Registered: 11-2010 |
Dear all, I hav a big problem in performing the S.Golay on my spectral data. The table is uploded to the Rcmdr, I also uploded the pakage RTisean, which has S.Golay,but the command that is written in the script is incomprehensible for me..please help! Magda |