[R] mb.long {timecourse}
How do you identify the genes which are differentially expressed using the mb.long function? More specifically, in the fruitfly example (see below), we begin with an expression matrix containing 2000 genes. How do I obtain not only the proportion of differentially expressed genes, but also the subset of genes which are differentially expressed? data(fruitfly) colnames(fruitfly) ## check if arrays are arranged in the default order gnames - rownames(fruitfly) assay - rep(c(A, B, C), each = 12) time.grp - rep(c(1:12), 3) size - rep(3, nrow(fruitfly)) out1 - mb.long(fruitfly, times=12, reps=size, rep.grp = assay, time.grp = time.grp) summary(out1) plotProfile(out1, type=b, gnames=gnames, legloc=c(2,15), pch=c(A,B,C), xlab=Hour) summary(out1) Length Class Mode M 72000 -none- numeric prop1 -none- numeric nu 1 -none- numeric Lambda1 121 -none- numeric percent 2 -none- numeric size 2000 -none- numeric con.group 36 -none- numeric rep.group 36 -none- numeric time.group 36 -none- numeric HotellingT2 2000 -none- numeric pos.HotellingT2 2000 -none- numeric Thanks to anyone who can shed light on this matter. -- View this message in context: http://r.789695.n4.nabble.com/mb-long-timecourse-tp4515669p4515669.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] limma design matrix
I am trying to construct my design matrix needed in the {limma} function lmfit but am having trouble with the formula I am to specify in the function model.matrix. Namely when to I use ~0 + factors (ex 1) vs ~-1 + factors (ex 2). Any clarification on this would be greatly appreciated. Thanks in advanced. Ex 1: f - factor(targets$Target, levels=c(RNA1,RNA2,RNA3)) design - model.matrix(~0+f) [page 409, of Bioinformatics and computational biology solutions using R and Bioconductor] and Ex 2: library(ALL) pdat - pData(ALL) design - model.matrix(~-1 + factor(pdat$type) [page 237, of Bioinformatics and computational biology solutions using R and Bioconductor] -- View this message in context: http://r.789695.n4.nabble.com/limma-design-matrix-tp4503651p4503651.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] hypergeometric function in ‘ mvtnorm’
Thanks for your advice. I actually meant to ask about the pmvt for the distribution function. Viewing the source code pmvt uses the function mvt which uses the function probval which sources the fortran code: Fortran(mvtdst, N = as.integer(n), NU = as.integer(df), LOWER = as.double(lower), UPPER = as.double(upper), INFIN = as.integer(infin), CORREL = as.double(corrF), DELTA = as.double(delta), MAXPTS = as.integer(x$maxpts), ABSEPS = as.double(x$abseps), RELEPS = as.double(x$releps), error = as.double(error), value = as.double(value), inform = as.integer(inform), PACKAGE = mvtnorm) I wish to look at how this mvtdst calculates the hypergeometric function (2_F_1). Anyway that I can see that? Thanks -- View this message in context: http://r.789695.n4.nabble.com/hypergeometric-function-in-mvtnorm-tp4483730p4485277.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] hypergeometric function in ‘ mvtnorm’
Is there any way to know how the dmvt function computes the hypergeometric function needed in the calculation for the density of multivariate t distribution? -- View this message in context: http://r.789695.n4.nabble.com/hypergeometric-function-in-mvtnorm-tp4483730p4483730.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] try to silence errors
I am trying to use the dmt function in the package {mnormt}. Throughout my algorithm, the covariance matrix is sometime calculated to be singular. When attempting to calculate the dmt function with a covariance that is not positive definite, I would like it to return Inf or NaN instead of an error message. I have been using the try function, however it is not yeilding the desired result. (I did not include the values since they are 17 dimensions, all you need to know is that sigij is singular) try(out - dmt(y[,j],new$mu[,,4],sigij,ceiling(new$nu[4])),silent=TRUE) Error in chol.default(x, pivot = FALSE) : the leading minor of order 16 is not positive definite I have similarly tried the trycatch function (see below for code) however I am using this within a function made my be (which i call ll) so I get errors like ll - function(new,y){ ... out - tryCatch(dmt(y[,j],new$mu[,,3],sigij, ceiling(new$nu[3])),error = function(e) -Inf) } ll(new,y) Error in chol.default(x, pivot = FALSE) : the leading minor of order 16 is not positive definite [1] -Inf If anyone could shed some light on this problem, it would be greatly appreciated. -- View this message in context: http://r.789695.n4.nabble.com/try-to-silence-errors-tp4225754p4225754.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] breaking up n object into random groups
say n = 100 I want to partition this into 4 random groups wheren n1 + n2 + n3 + n4 = n and ni is the number of elements in group i. Thank you for you help -- View this message in context: http://r.789695.n4.nabble.com/breaking-up-n-object-into-random-groups-tp4147476p4147476.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] tryCatch with integration
I am trying to compute the approximate numerical integration of the following expression using the integrate function: (integrate(function(x) {log(1+x^2)*(1+x^2)^(-20.04543)},low,Inf)$val) Error in integrate(function(x) { : the integral is probably divergent which gives me an error. If this error happens I want to return the following instead (area(function(x) {-2*log(cos(x))*(cos(x)^(2*e-2))}, atan(low), pi/2)) [1] 0.01045636 In order to accomplish this I tried the following: this= tryCatch((integrate(function(x) {log(1+x^2)*(1+x^2)^(-e)},low,Inf)$val),error=function(e) e, finally = (area(function(x) {-2*log(cos(x))*(cos(x)^(2*e-2))}, atan(low), pi/2))) this simpleError in integrate(function(x) {log(1 + x^2) * (1 + x^2)^(-e)}, low, Inf): the integral is probably divergent Perhaps I am using tryCatch incorrectly? Any help would be greatly appreciated. -- View this message in context: http://r.789695.n4.nabble.com/tryCatch-with-integration-tp3954528p3954528.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] symmetric matrix multiplication
I have a symmetric matrix B (17x17), and a (17x17) square matrix A. If do the following matrix multiplication I SHOULD get a symmetric matrix, however i don't. The computation required is: C = t(A)%*%B%*%A here are some checks for symmetry (max(abs(B - t(B [1] 0 C = t(A)%*%B%*%A (max(abs(C - t(C [1] 3.552714e-15 Any help on the matter would be very much appreciated. -- View this message in context: http://r.789695.n4.nabble.com/symmetric-matrix-multiplication-tp3929739p3929739.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] symmetric matrix multiplication
Thank you Dan and Ted for these helpful comments. I will implement this simple force symmetry code you suggested and make sure I familiarize with this floating-point calculation problem so I can recognize such issues in the future. -- View this message in context: http://r.789695.n4.nabble.com/symmetric-matrix-multiplication-tp3929739p3930832.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] sampling from the multivariate truncated normal
Well, for 0.828324 x[2] Inf the probablility is roughly 0 hence not easy to draw random numbers out there Uwe Ligges How is this probability roughly 0? -- View this message in context: http://r.789695.n4.nabble.com/sampling-from-the-multivariate-truncated-normal-tp3626438p3647039.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] sampling from the multivariate truncated normal
I am trying generate a sample for a truncated multivariate normal distribution via the rtmvnorm function in the {tmvtnorm} package. Why does the following produce NaNs? rtmvnorm(1, mean = rep(0, 2), matrix(c(0.06906084, -0.07463565, -0.07463565, 0.08078086),2),c(-0.4316738, 0.8283240), c(Inf,Inf), algorithm=gibbsR, burn.in.samples=100) Thanks -- View this message in context: http://r.789695.n4.nabble.com/sampling-from-the-multivariate-truncated-normal-tp3626438p3626438.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] rtmvt
I want to use the rtmvt from the {tmvtnorm} package using the gibbs algorithm but how to i specify the nested function rtmvnorm to use gibbs as well? Right now I am using the code: for (i in 1:g){ for (j in 1:n){ sgamma[,,i,j] = rtmvt(n=50, mean=mu[i,j], sigma[i,j], df=nu[i], lower=rep(0,2),algorithm=gibbs) } } heres an example of one iteration: mu[1,1] -0.09734357 0.51578628 sigma[1,1] [,1] [,2] [1,] 0.4250681 0.0253649 [2,] 0.0253649 0.4250681 when I run this i get 50 errors saying: Warning messages: 1: In rtmvnorm.rejection(n, mean, sigma, lower, upper, ...) : Acceptance rate is very low and rejection sampling becomes inefficient. Consider using Gibbs sampling. I have figured out that this is coming from the internal function rtmvnorm who's default is the rejection algorithm. Is there any way that I can specify that to be Gibbs as well? Thanks -- View this message in context: http://r.789695.n4.nabble.com/rtmvt-tp3562217p3562217.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] rtmvt
I have been using the rtmvt function in the {tmvtnorm} package i'm getting the warning: Acceptance rate is very low and rejection sampling becomes inefficient. Consider using Gibbs sampling. but i AM specifying the gibbs algorithm!!: rtmvt(M, mean=q[,,i,j], sigma=((u[i,j] + nu[i])/(p+nu[i]))*delta[,,i], df=ceiling(nu[i]+p), lower=c(0,0), algorithm=gibbs) Any ideas why I am getting this warning and how I can fix it? -- View this message in context: http://r.789695.n4.nabble.com/rtmvt-tp3440751p3440751.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] multivariate t distribution
I have been working the the pmt function in the {mnormt} package and which requires S a positive definite matrix representing the scale matrix of the distribution, such that S*df/(df-2) is the variance-covariance matrix when df2; a vector of length 1 is also allowed (in this case, d=1 is set) is there a way that I can specify the scale covariance matrix instead? Or alternatively, how do I convert the scale covariance matrix into this positive definite S matrix. Thanks in advanced. -- View this message in context: http://r.789695.n4.nabble.com/multivariate-t-distribution-tp3432257p3432257.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] truncated distributions
I am sampling from the truncated multivariate student t distribution rtmvt in the package {tmvtnorm}. My question is about the mean vector. Is it possible to define a mean vector outside of the truncated region? Thank you in advance for any help. -- View this message in context: http://r.789695.n4.nabble.com/truncated-distributions-tp3422245p3422245.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] truncated distributions
The definition of the mean vector is essentially what my question boils down to. In the functions details, the author states We sample x ~ T(mean, Sigma, df) subject to the rectangular truncation lower = x = upper. Currently, two random number generation methods are implemented: rejection sampling and the Gibbs Sampler. So if the mean vector in the rtmvt function is the mean of the parent distribution's mean (as I hope it is), then it would be acceptable to define a mean vector outside of the truncated range. Clarification of this point would be greatly appreciated. -- View this message in context: http://r.789695.n4.nabble.com/truncated-distributions-tp3422245p3422434.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] pmt
I am working with the pmt function in the {mnormt} package, and i am getting negative values returned. the following is an example of one of my outputs: pmt(x = c(3.024960, -1.010898), mean = c(21.18844, 21.18844), S = matrix(c(.319,.139,.139,0.319), 2, 2),df = 42) # -6.585641e-18 Any help on why i'm getting negative numbers would be very much appreciated. THanks! -- View this message in context: http://r.789695.n4.nabble.com/pmt-tp3409543p3409543.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.