Re: [R] lmom package - Resending the email
Katherine, for a deeper understanding of differing values it makes sense to provide the list at least with an online description of the corresponding functions used in Minitab and SPSS… Best Simon On 03 Dec 2014, at 10:45, Katherine Gobin via R-help r-help@r-project.org wrote: Dear R forum I sincerely apologize as my earlier mail with the captioned subject, since all the values got mixed up and the email is not readable. I am trying to write it again. My problem is I have a set of data and I am trying to fit some distributions to it. As a part of this exercise, I need to find out the parameter values of various distributions e.g. Normal distribution, Log normal distribution etc. I am using lmom package to do the same, however the parameter values obtained using lmom pacakge differ to a large extent from the parameter values obtained using say MINITAB and SPSS as given below - _ amounts = c(38572.5599129508,11426.6705314315,21974.1571641187,118530.32782443,3735.43055996748,66309.5211176106,72039.2934132668,21934.8841708626,78564.9136114375,1703.65825161293,2116.89180930203,11003.495671332,19486.3296339113,1871.35861218795,6887.53851253407,148900.978055447,7078.56497101651,79348.1239806592,20157.6241066905,1259.99802108593,3934.45912233674,3297.69946631591,56221.1154121067,13322.0705174134,45110.2498756567,31910.3686613912,3196.71168501252,32843.0140437202,14615.1499458453,13013.9915051561,116104.176753387,7229.03056392023,9833.37962177814,2882.63239493673,165457.372543821,41114.066453219,47188.1677766245,25708.5883755617,82703.7378298092,8845.04197017415,844.28834047836,35410.8486123933,19446.3808445684,17662.2398792892,11882.8497070776,4277181.17817307,30239.0371267968,45165.7512343364,22102.8513746687,5988.69296597127,51345.0146170238,1275658.35495898,15260.4892854214,8861.76578480635,37647.1638704867,4979.53544046949,7012.48134772332,3385.20612391205,1911.03114395959,66886.5036605189,2223.47536156462,814.947809578378,234.028589468841,5397.4347625133,13346.3226579065,28809.3901352898,6387.69226236731,5639.42730553242,2011100.92675507,4150.63707173462,34098.7514446498,3437.10672573502,289710.315303182,8664.66947305203,13813.3867161134,208817.521491857,169317.624400274,9966.78447705792,37811.1721605562,2263.19211279927,80434.5581206454,19057.8093104899,24664.5067589624,25136.5042354789,3582.85741610706,6683.13898432794,65423.9991390846,134848.302304064,3018.55371579808,546249.641168158,172926.689143006,3074.15064180208,1521.70624812788,59012.4248281661,21226.928522236,17572.5682970983,226.646947337851,56232.2982652019,14641.0043361533,6997.94414914865) library(lmom) lmom = samlmu(amounts) # __ # Normal Distribution parameters parameters_of_NOR - pelnor(lmom); parameters_of_NOR mu sigma 115148.4175945.8 Location Scale Minitab 115148.4 485173SPSS 115148.4 485173 # __ # Log Normal (3 Parameter) Distribution parameters zetamu sigma 3225.7988909.114879 2.240841 LocationScale Shape MINITAB 9.73361 1.76298 75.51864SPSS 9.73361.763 75.519 # __ Besides Genaralized extreme Value distributions, all the other distributions e.g. Gamma, Exponential (2 parameter) distributions etc give different results than MINITAB and SPSS. Can some one guide me? Regards Katherine [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] Building R for better performance
calculation (vector calc) 1.310.38 Creation of a 6000x6000 Hilbert matrix (matrix calc) 0.77 0.99 Grand common divisors of 400,000 pairs (recursion) 0.63 0.56 Creation of a 1000x1000 Toeplitz matrix (loops) 2.24 2.34 Escoufier's method on a 90x90 matrix (mixed) 9.55 6.02 Total 274.93 21.01 Regards, Jonathan Anspach Sr. Software Engineer Intel Corp. jonathan.p.ansp...@intel.com 713-751-9460 -Original Message- From: Simon Zehnder [mailto:szehn...@uni-bonn.de] Sent: Wednesday, March 05, 2014 3:55 AM To: Anspach, Jonathan P Cc: r-help@r-project.org Subject: Re: [R] Building R for better performance Jonathan, I myself tried something like this - comparing gcc, clang and intel on a Mac. From my experiences in HPC on the university cluster (where we also use the Xeon Phi, Landeshochleistungscluster University RWTH Aachen), the Intel compiler has better code optimization in regard to vectorisation, etc. (clang is up to now suffering from a not yet implemented OpenMP library). Here is a revolutionanalytics article about this topic: http://blog.revolutionanalytics.com/2010/06/performance-benefits-of-multithreaded-r.html As I usually use the Rcpp package for C++ extensions this could give me further performance. Though, I already failed when trying to compile R with the Intel compiler and linking against the MKL (see my topic in the Intel developer zone: http://software.intel.com/en-us/comment/1767418 and my threads on the R-User list: https://stat.ethz.ch/pipermail/r-sig-mac/2013-November/010472.html). So, to your questions: 1) I think that most admins do not even use the Intel compiler to compile R - this seems to me rare. There are some people I know they do and I think they could be aware of it - but these are only a few. As R is growing in usage and I do know from regional user meetings that very large companies start using it in their BI units - this should be of interest. 2) I would really welcome this step because compilation with intel (especially on a Mac) and linking to the MKL seems to be delicate. I am interested in the data - so if it is possible send it via the list or directly to my account. Further, could you show some code that you used for the computations? Best Simon On 04 Mar 2014, at 22:44, Anspach, Jonathan P jonathan.p.ansp...@intel.com wrote: Greetings, I'm a software engineer with Intel. Recently I've been investigating R performance on Intel Xeon and Xeon Phi processors and RH Linux. I've also compared the performance of R built with the Intel compilers and Intel Math Kernel Library to a default build (no config options) that uses the GNU compilers. To my dismay, I've found that the GNU build always runs on a single CPU core, even during matrix operations. The Intel build runs matrix operations on multiple cores, so it is much faster on those operations. Running the benchmark-2.5 on a 24 core Xeon system, the Intel build is 13x faster than the GNU build (21 seconds vs 275 seconds). Unfortunately, this advantage is not documented anywhere that I can see. Building with the Intel tools is very easy. Assuming the tools are installed in /opt/intel/composerxe, the process is simply (in bash shell): $ . /opt/intel/composerxe/bin/compilervars.sh intel64 $ ./configure --with-blas=-L/opt/intel/composerxe/mkl/lib/intel64 -lmkl_intel_lp64 -lmkl_intel_thread -lmkl_core -liomp5 -lpthread -lm --with-lapack CC=icc CFLAGS=-O2 CXX=icpc CXXFLAGS=-O2 F77=ifort FFLAGS=-O2 FC=ifort FCFLAGS=-O2 $ make $ make check My questions are: 1) Do most system admins and/or R installers know about this performance difference, and use the Intel tools to build R? 2) Can we add information on the advantage of building with the Intel tools, and how to do it, to the installation instructions and FAQ? I can post my data if anyone is interested. Thanks, Jonathan Anspach Sr. Software Engineer Intel Corp. jonathan.p.ansp...@intel.com 713-751-9460 __ 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-benchmark-25_large.R
Re: [R] Building R for better performance
Jonathan, I myself tried something like this - comparing gcc, clang and intel on a Mac. From my experiences in HPC on the university cluster (where we also use the Xeon Phi, Landeshochleistungscluster University RWTH Aachen), the Intel compiler has better code optimization in regard to vectorisation, etc. (clang is up to now suffering from a not yet implemented OpenMP library). Here is a revolutionanalytics article about this topic: http://blog.revolutionanalytics.com/2010/06/performance-benefits-of-multithreaded-r.html As I usually use the Rcpp package for C++ extensions this could give me further performance. Though, I already failed when trying to compile R with the Intel compiler and linking against the MKL (see my topic in the Intel developer zone: http://software.intel.com/en-us/comment/1767418 and my threads on the R-User list: https://stat.ethz.ch/pipermail/r-sig-mac/2013-November/010472.html). So, to your questions: 1) I think that most admins do not even use the Intel compiler to compile R - this seems to me rare. There are some people I know they do and I think they could be aware of it - but these are only a few. As R is growing in usage and I do know from regional user meetings that very large companies start using it in their BI units - this should be of interest. 2) I would really welcome this step because compilation with intel (especially on a Mac) and linking to the MKL seems to be delicate. I am interested in the data - so if it is possible send it via the list or directly to my account. Further, could you show some code that you used for the computations? Best Simon On 04 Mar 2014, at 22:44, Anspach, Jonathan P jonathan.p.ansp...@intel.com wrote: Greetings, I'm a software engineer with Intel. Recently I've been investigating R performance on Intel Xeon and Xeon Phi processors and RH Linux. I've also compared the performance of R built with the Intel compilers and Intel Math Kernel Library to a default build (no config options) that uses the GNU compilers. To my dismay, I've found that the GNU build always runs on a single CPU core, even during matrix operations. The Intel build runs matrix operations on multiple cores, so it is much faster on those operations. Running the benchmark-2.5 on a 24 core Xeon system, the Intel build is 13x faster than the GNU build (21 seconds vs 275 seconds). Unfortunately, this advantage is not documented anywhere that I can see. Building with the Intel tools is very easy. Assuming the tools are installed in /opt/intel/composerxe, the process is simply (in bash shell): $ . /opt/intel/composerxe/bin/compilervars.sh intel64 $ ./configure --with-blas=-L/opt/intel/composerxe/mkl/lib/intel64 -lmkl_intel_lp64 -lmkl_intel_thread -lmkl_core -liomp5 -lpthread -lm --with-lapack CC=icc CFLAGS=-O2 CXX=icpc CXXFLAGS=-O2 F77=ifort FFLAGS=-O2 FC=ifort FCFLAGS=-O2 $ make $ make check My questions are: 1) Do most system admins and/or R installers know about this performance difference, and use the Intel tools to build R? 2) Can we add information on the advantage of building with the Intel tools, and how to do it, to the installation instructions and FAQ? I can post my data if anyone is interested. Thanks, Jonathan Anspach Sr. Software Engineer Intel Corp. jonathan.p.ansp...@intel.com 713-751-9460 __ 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-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] Where did lost variables go
A reproducible example would do well here David Best Simon On 31 Dec 2013, at 02:42, David Parkhurst parkh...@indiana.edu wrote: I have several variables in a data frame that aren't listed by ls() after I attach that data frame. Where did they go, and how can I stop the hidden ones from masking the local ones? Thanks for any help. David [[alternative HTML version deleted]] __ 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-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] Estimation of AR(1) model by QML.
Why not using optim on the likelihood in a) with normally distributed standard errors and for b) optim with a likelihood with t(3)-distributed standard errors? Best Simon On 30 Dec 2013, at 21:19, Xuse Chuse chus...@gmail.com wrote: Dear Users, I am trying to estimate a model Y(t)=alpha+rho*Y(t-1)+e(t) where i know e(t)~t(3). a) I want to estimate (alpha, rho) by QML estimation assuming (wrongly) that e(t)~N(0,sigma2) and calculate the standard errors. b) Estimate (alpha, rho) by ML estimation assuming (correctly) e(t)~t(3) and compute standard errors. Can anyone help me out to figure out how I could answer this question in R? Thank you beforehand. Cheers, Chuse [[alternative HTML version deleted]] __ 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-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] Milliken and Johnson Unbalanced Machines
I don’t know anything about this topic, but looking into this book: http://books.google.de/books?id=3TVDQBAJpg=PA25lpg=PA25dq=milliken+and+johnson+unbalanced+machinessource=blots=lxqStiQju5sig=d5dG_cHsTzilCIklBrW9SAMKYRMhl=desa=Xei=ifrBUuvoM8qPtAbPqoCQDQved=0CGUQ6AEwBQ#v=onepageq=milliken%20and%20johnson%20unbalanced%20machinesf=false should give you a little hint I guess. Best Simon On 30 Dec 2013, at 23:37, Troels Ring tr...@gvdnet.dk wrote: Dear friends - reading Milliken and Johnson on messy data I failed to find R code to master the unbalanced Machines data in ch 23. I wonder if anyone can lend me a hand - again Happy new year Troels Ring Aalborg, Denmark __ 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-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] testing for bimodal and for dip in between modes in R
Testing for bimodality is rather testing for unimodality. Hartigan and Hartigan (1985) presented the Dip-Test which is implemented in the R package DipTest with a much better approximation of the test distribution. If the test statistic is too high unimodality is rejected. To estimate the dip point you could choose among several possibilities: (1) A very easy method is to use the kmeans function for a kmeans cluster and use the point in the middle of the connecting line between the kmeans cluster centers. (2) You could estimate a finite mixture distribution and take the middle of the connecting line of the modes. Best Simon On 24 Nov 2013, at 20:41, Felix Breden bre...@sfu.ca wrote: Hi I have distributions that are typically bimodal (see attached .pdf), and I would like to test for bimodality, and then estimate the point between the two modes, the dip in the distributions. any help would be greatly appreciated. thanks felix m66.junction.aln.pairwise.histogram.pdf__ 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-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] MM estimator
Thanks Martin, for making this clear to me, I thought of Pearson’s Method-of-Moments. On 14 Nov 2013, at 16:08, Martin Maechler maech...@stat.math.ethz.ch wrote: SZ == Simon Zehnder szehn...@uni-bonn.de on Thu, 14 Nov 2013 08:52:16 +0100 writes: SZ Check the gmm package with a weighting matrix equal to SZ the identity. Best SZ Simon No.MM is ambigous: The gmm package is about Generalized Method of Moments but Izhak is really asking about MM estimation in the field of robust statistics, where an MM estimator is a special kind of M-Estimator (and these, introduced by Peter Huber (1964), Annals, fortunately are well known unambigously, as a generalization of ML estimators (= MLE)), namely an M-estimator with redescending psi (== non-convex problem == solution typically depends on starting value), started with a high-breakdown point initial estimate. Izhak mentioned nlrob() from package robustbase, and he is right that this does not provide an MM estimator, but only an M estimate started by Least Squares. Fortunately, we have had contributions to the robustbase package from Eduardo Conceição since last summer. He provided alternative versions of nlrob(), notably a nlrob.MM() currently *hidden* {and undocumented}. I.e. you need install.packages(robustbase, repos=http://R-Forge.R-project.org;) library(robustbase) and then use robustbase:::nlrob.MM() The reason for all this is that I plan to have one nlrob() function but with a new argument 'method' and so eventually nlrob.MM(, method=MM) will correspond to current to the R-forge version robustbase:::nlrob.MM() --- Please for all more, use the dedicated mailing list R-SIG-robust only (i.e. do *not* just reply-all to this e-mail!). Martin Maechler, ETH Zurich SZ On 14 Nov 2013, at 04:37, IZHAK shabsogh SZ ishaqb...@yahoo.com wrote: hi I have a nonlinear regression model with two parameter, i have estimated using nls in R and i want to also find the estimate using MM, someone refer me to this function nlrob but this is main for only M estimate. can you please help me findout a function in R for MM nonlinear regression thanks [[alternative HTML version deleted]] __ 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. SZ __ SZ R-help@r-project.org mailing list SZ https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do SZ read the posting guide SZ http://www.R-project.org/posting-guide.html and provide SZ commented, minimal, self-contained, reproducible code. __ 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] optimization: multiple assignment problem
It would be more clear if you tell, what you want to do instead of what you do not want to do. If you start with a usual cost matrix (whatever cost function you have) and you have to assign N to N this reduces to the well-known Munkre’s algorithm (see for example: http://gallery.rcpp.org/articles/minimal-assignment/). If you have to assign M to N, it is an extended version of the Munkre’s algorithm which has been covered by Bourgeois and Lassalle (1971). All this is graph theory. Best Simon On 14 Nov 2013, at 19:24, Jean-Francois Chevalier jean-francois.cheval...@bisnode.com wrote: Hello, I'm trying to solve a multiple assignment problem. I found a package Adagio and its function mknapsack which maximize vstar = p(1)*(x(1,1) + ... + x(m,1)) + ... ... + p(n)*(x(1,n) + ... + x(m,n)) subject to w(1)*x(i,1) + ... + w(n)*x(i,n) = k(i) for i=1,...,m x(1,j) + ... + x(m,j) = 1 for j=1,...,n x(i,j) = 0 or 1 for i=1,...,m , j=1,...,n , It's close to what I'm trying to do except that 1)k(i) = k for any I (not an issue) 2)p is dependent of the item AND the knapsack 3)each item must be assigned maximize vstar = p(1,1)*x(1,1) + ... + p(m,1)*x(m,1) + ... ... + p(1,n)*x(1,n) + ... + p(m,n)*x(m,n) with p(j,i) profit of assigning item i to knapsack j subject to w(1)*x(i,1) + ... + w(n)*x(i,n) = k for i=1,...,m x(1,j) + ... + x(m,j) = 1 for j=1,...,n x(i,j) = 0 or 1 for i=1,...,m , j=1,...,n , It would be really helpful if you could indicate me any package, function that would solve my problem? Thanks in advance, Best regards, Jean-François DISCLAIMER \ This e-mail and any attachments t...{{dropped:13}} __ 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-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] issues with calling predict.coxph.penal (survival) inside a function - subset-vector not found. Because of NextMethod?
Works for me: predicting_function(fit2,test1) 1 2 3 4 5 6 7 -1.0481141 0.1495946 0.4492597 0.4492597 0.9982492 -0.4991246 -0.4991246 Best Simon On 13 Nov 2013, at 15:46, julian.bo...@elitepartner.de wrote: Hello everyone, I got an issue with calling predict.coxph.penal inside a function. Regarding the context: My original problem is that I wrote a function that uses predict.coxph and survfit(model) to predict a lot of survival-curves using only the basis-curves for the strata (as delivered by survfit(model) ) and then adapts them with the predicted risk-scores. Because there are cases where my new data has strata which didn't exist in the original model I exclude them, using a Boolean vector inside the function. I end up with a call like this: predict (coxph_model, newdata[subscript_vector,] ) This works fine for coxph.model, but when I fit a model with a spline (class coxph.penal), I get an error: Error in `[.data.frame`(newdata, [subscript_vector, ) : object '[subscript_vector ' not found I suppose this is because of NextMethod, but I am not sure how to work around it. I also read a little bit about all those matching-and-frame-issues, But must confess I am not really into it. I attach a reproducible example. Any help or suggestions of work-arounds will be appreciated. Thanks Julian version _ platform x86_64-w64-mingw32 arch x86_64 os mingw32 system x86_64, mingw32 status major 3 minor 0.1 year 2013 month 05 day16 svn rev62743 language R version.string R version 3.0.1 (2013-05-16) nickname Good Sport ##TEST-DATA # Create the simplest test data set test1 - data.frame(time=c(4,3,1,1,2,2,3), status=c(1,1,1,0,1,1,0), x=c(0,2,1,1,1,0,0), sex=c(0,0,0,0,1,1,1)) # Fit a stratified model fit1 - coxph(Surv(time, status) ~ x + strata(sex), test1) summary(fit1) #fit stratified wih spline fit2 - coxph(Surv(time, status) ~ pspline(x, df=2) + strata(sex), test1) summary(fit2) #function to predict within predicting_function - function(model, newdata){ subs -vector(mode='logical', length=nrow(newdata)) subs[1:length(subs)]- TRUE ret - predict (model, newdata=newdata[subs,]) return(ret) } predicting_function(fit1, test1) # works predicting_function(fit2,test1) #doesnt work - Error in `[.data.frame`(newdata, subs, ) : object 'subs' not found # probably because of NextMethod # traceback() #12: `[.data.frame`(newdata, subs, ) #11: newdata[subs, ] #10: is.data.frame(data) #9: model.frame.default(data = newdata[subs, ], formula = ~pspline(x, # df = 2) + strata(sex), na.action = function (object, ...) # object) #8: model.frame(data = newdata[subs, ], formula = ~pspline(x, df = 2) + # strata(sex), na.action = function (object, ...) # object) #7: eval(expr, envir, enclos) #6: eval(tcall, parent.frame()) #5: predict.coxph(model, newdata = newdata[subs, ]) #4: NextMethod(predict, object, ...) #3: predict.coxph.penal(model, newdata = newdata[subs, ]) #2: predict(model, newdata = newdata[subs, ]) at #5 #1: predicting_function(fit2, test1) [[alternative HTML version deleted]] __ 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-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] MM estmator
Check the gmm package with a weighting matrix equal to the identity. Best Simon On 14 Nov 2013, at 04:37, IZHAK shabsogh ishaqb...@yahoo.com wrote: hi I have a nonlinear regression model with two parameter, i have estimated using nls in R and i want to also find the estimate using MM, someone refer me to this function nlrob but this is main for only M estimate. can you please help me findout a function in R for MM nonlinear regression thanks [[alternative HTML version deleted]] __ 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-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] S4; Setter function is not chaning slot value as expected
If you want to set a slot you have to refer to it: ac@CustomerID - “54321” or you use your setter: ac - CustomerID(ac, “54321”) What you did was creating a new symbol CustomerID referring to the String “54321” CustomerID - “54321” CustomerID [1] “54321” Best Simon On 09 Nov 2013, at 15:31, daniel schnaider dschnai...@gmail.com wrote: It is my first time programming with S4 and I can't get the setter fuction to actually change the value of the slot created by the constructor. I guess it has to do with local copy, global copy, etc. of the variable - but, I could't find anything relevant in documentation. Tried to copy examples from the internet, but they had the same problem. # The code setClass (Account , representation ( customer_id = character, transactions = matrix) ) Account - function(id, t) { new(Account, customer_id = id, transactions = t) } setGeneric (CustomerID-, function(obj, id){standardGeneric(CustomerID-)}) setReplaceMethod(CustomerID, Account, function(obj, id){ obj@customer_id - id obj }) ac - Account(12345, matrix(c(1,2,3,4,5,6), ncol=2)) ac CustomerID - 54321 ac #Output ac An object of class Account Slot customer_id: [1] 12345 Slot transactions: [,1] [,2] [1,]14 [2,]25 [3,]36 # CustomerID is value has changed to 54321, but as you can see it does't CustomerID - 54321 ac An object of class Account Slot customer_id: [1] 12345 Slot transactions: [,1] [,2] [1,]14 [2,]25 [3,]36 Help! [[alternative HTML version deleted]] __ 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-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] S4 vs S3. New Package
This depends very often of on the developer and what he is comfortable with. I like S4 classes, as I come from C++ and S4 classes approximate C++ classes at least more than S3 classes do (Reference Classes would do so even more and I know very good R programmers liking these most). 1) I wrote a package for MCMC simulation with S4 classes carrying all simulated values - fast enough for me: in less than 1.5 secs I have my sample of 100.000 values together with several other 100T values like log-likelihoods, posterior hyper parameters, etc. I watch out for not copying too often an object but sometimes it is not avoidable. 2) That is not true: Books: http://www.amazon.de/Software-Data-Analysis-Programming-Statistics/dp/0387759352/ref=sr_1_1?ie=UTF8qid=1384014486sr=8-1keywords=John+chambers+data http://www.amazon.de/Programming-Data-Language-John-Chambers/dp/0387985034/ref=sr_1_4?ie=UTF8qid=1384014486sr=8-4keywords=John+chambers+data Online: https://www.rmetrics.org/files/Meielisalp2009/Presentations/Chalabi1.pdf https://www.stat.auckland.ac.nz/S-Workshop/Gentleman/S4Objects.pdf And for a bunch of packages look into the Bioconductor packages. Best Simon On 09 Nov 2013, at 16:22, daniel schnaider dschnai...@gmail.com wrote: Hi, I am working on a new credit portfolio optimization package. My question is if it is more recommended to develop in S4 object oriented or S3. It would be more naturally to develop in object oriented paradigm, but there is many concerns regarding S4. 1) Performance of S4 could be an issue as a setter function, actually changes the whole object behind the scenes. 2) Documentation. It has been really hard to find examples in S4. Most books and articles consider straightforward S3 examples. Thanks, [[alternative HTML version deleted]] __ 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-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] select .txt from .txt in a directory
I do not understand the question. If you already know the names what is the problem to select the files by names? If you have the names but not inside of R you have to find a name pattern to avoid typing them in. Is there a pattern, e.g. da.txt, db.txt, dc.txt? On 08 Nov 2013, at 18:33, Zilefac Elvis zilefacel...@yahoo.com wrote: Hi, I have 300 .txt files in a directory. Out of this 300, I need just 100 of the files. I have the names of the 100 .txt files which are also found in the 300 .txt files. How can I extract only the 100 .txt files from the 300 ,txt files? e.g given d1.txt, ds.txt, dx.txt, df.txt...d300.txt, how can I select only d1.txt and df.txt? Remember, I have 300 of such and want to extract 100 of them with names known. Thanks for your great help. Atem. [[alternative HTML version deleted]] __ 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-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] select .txt from .txt in a directory
If you want to type in the names by hand, you can simply use read.table to load them into R … I still don’t get the aim of your text file handling On 08 Nov 2013, at 18:51, Zilefac Elvis zilefacel...@yahoo.com wrote: All files are text files. They are found in a folder on my computer. Assume that I know the names of some of the files I want to select from the 300 txt files. How can I do this in R. Atem. On Friday, November 8, 2013 11:44 AM, Simon Zehnder szehn...@uni-bonn.de wrote: I do not understand the question. If you already know the names what is the problem to select the files by names? If you have the names but not inside of R you have to find a name pattern to avoid typing them in. Is there a pattern, e.g. da.txt, db.txt, dc.txt? On 08 Nov 2013, at 18:33, Zilefac Elvis zilefacel...@yahoo.com wrote: Hi, I have 300 .txt files in a directory. Out of this 300, I need just 100 of the files. I have the names of the 100 .txt files which are also found in the 300 .txt files. How can I extract only the 100 .txt files from the 300 ,txt files? e.g given d1.txt, ds.txt, dx.txt, df.txt...d300.txt, how can I select only d1.txt and df.txt? Remember, I have 300 of such and want to extract 100 of them with names known. Thanks for your great help. Atem. [[alternative HTML version deleted]] __ 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-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] select .txt from .txt in a directory
Elvis, first, keep things on the list - so others can learn and comment. Second, as Sarah already commented: We do not like to open unsolicited binary attachments on the list. Sarah gives a good hint how to post data to the list. What I would do so far is use the matching columns to get the names you need from BTemperature: temp_inv - read.table(Temperature Inventory, … ) (here I would change the .xlsx to a .csv and use read.csv instead of read.table) btemp - read.table(“BTemperature_Stations.txt”, … ) (again think about converting via Excel to .csv - it makes things far more easy) Check ?read.table for options - you gonna need them. Then match mynames - btemp[(temp_inv[,3] %in% btemp[, 3]), 2] Now you have the names of the stations and if your .txt files are named by the stations you can do something like: for (name in mynames) { tmp.table - read.table(paste(“path/to/your/Homog_daily_min_temp/“, name, “.txt”, sep = “”), … ) …. do things with the data } Best Simon On 08 Nov 2013, at 19:26, Zilefac Elvis zilefacel...@yahoo.com wrote: Hi Simon, Attached are my data files. Btemperature_Stations is my main file. Temperature inventory is my 'wanted' file and is a subset of Btemperature_Stations. Using column 3 in both files, select the files in Temperature inventory from Btemperature_Stations. The .zip file contains the .txt files which you will extract to a folder and do the selection in R. Thanks, Atem. On Friday, November 8, 2013 11:54 AM, Simon Zehnder szehn...@uni-bonn.de wrote: If you want to type in the names by hand, you can simply use read.table to load them into R … I still don’t get the aim of your text file handling On 08 Nov 2013, at 18:51, Zilefac Elvis zilefacel...@yahoo.com wrote: All files are text files. They are found in a folder on my computer. Assume that I know the names of some of the files I want to select from the 300 txt files. How can I do this in R. Atem. On Friday, November 8, 2013 11:44 AM, Simon Zehnder szehn...@uni-bonn.de wrote: I do not understand the question. If you already know the names what is the problem to select the files by names? If you have the names but not inside of R you have to find a name pattern to avoid typing them in. Is there a pattern, e.g. da.txt, db.txt, dc.txt? On 08 Nov 2013, at 18:33, Zilefac Elvis zilefacel...@yahoo.com wrote: Hi, I have 300 .txt files in a directory. Out of this 300, I need just 100 of the files. I have the names of the 100 .txt files which are also found in the 300 .txt files. How can I extract only the 100 .txt files from the 300 ,txt files? e.g given d1.txt, ds.txt, dx.txt, df.txt...d300.txt, how can I select only d1.txt and df.txt? Remember, I have 300 of such and want to extract 100 of them with names known. Thanks for your great help. Atem. [[alternative HTML version deleted]] __ 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. BTemperature_Stations.txtTempearture inventory.xlsxHomog_daily_min_temp.zip __ 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] R with openblas and atlas
There is no R code following On 01 Nov 2013, at 05:34, Li Bowen bowenl...@gmail.com wrote: Hi, I have been trying to build R with optimized BLAS library. I am using a Ubuntu 13.10 x86_64 desktop, on which I am able to build R with openblas without any problem: #BEGIN_SRC sh ./configure --enable-BLAS-shlib --enable-R-shlib LIBnn=lib --disable-nls --with-blas=-L/usr/lib/openblas-base/ -lopenblas --enable-memory-profiling make sudo make install #END_SRC However, on redhat 5.9, I am not able to install openblas. Firstly, there is no pre-built package, even for later version of redhat. Secondly, I am not able to build openblas from source, actually not even able to install newer gcc from source. I then tried to install ATLAS on redhat and build R with it. -t 2: 2 threads #BEGIN_SRC sh ../configure --shared -t 2 -b 64 -D c -DPentiumCPS=1600 --with-netlib-lapack-tarfile=/path-to-lapack-3.4.2.tgz make build make check make ptcheck make time make install # R ../configure --enable-BLAS-shlib --enable-R-shlib --disable-nls --enable-memory-profiling --with-blas=-L/usr/local/atlas/lib -lptf77blas -lpthread -latlas make sudo make install #END_SRC Installation is successful. However when I run the following code, only one thread is used. I have looked through lots of manuals and forums and couldn't find an answer. Please advise. Thanks a lot. -- Sincerely, Bowen __ 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-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] Combinations of values in two columns
You could use the data.table package require(data.table) DT - data.table(Friend1 = sample(LETTERS, 10, replace = TRUE), Friend2 = sample(LETTERS, 10, replace = TRUE), Indicator = 1) ALL - data.table(unique(expand.grid(DT))) setkey(ALL) OTHERS - ALL[!DT] OTHERS[, Indicator := 0] RESULT - rbind(DT, ALL) Best Simon On 01 Nov 2013, at 10:32, Thomas thomas.ches...@nottingham.ac.uk wrote: I have data that looks like this: Friend1, Friend2 A, B A, C B, A C, D And I'd like to generate some more rows and another column. In the new column I'd like to add a 1 beside all the existing rows. That bit's easy enough. Then I'd like to add rows for all the possible directed combinations of rows not included in the existing data. So for the above I think that would be: A, D D, A B, C C, B B, D C, A D, B D, C and then put a 0 in the column beside these. Can anyone suggest how to do this? I'm using R version 2.15.3. Thank you, Thomas Chesney This message and any attachment are intended solely for the addressee and may contain confidential information. If you have received this message in error, please send it back to me, and immediately delete it. Please do not use, copy or disclose the information contained in this message or in any attachment. Any views or opinions expressed by the author of this email do not necessarily reflect the views of the University of Nottingham. This message has been checked for viruses but the contents of an attachment may still contain software viruses which could damage your computer system, you are advised to perform your own checks. Email communications with the University of Nottingham may be monitored as permitted by UK legislation. __ 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-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] set difference between two data frames
You could e.g. take the data.table package (every data.table is a data.frame) and make a join: dt.x - data.table(x) dt.y - data.table(y) merge.xy - x[y, nomatch = 0] diff.xy - x[!merge.xy] On 31 Oct 2013, at 21:41, Yasin Gocgun entropy...@gmail.com wrote: Thanks. Actually, I forgot to add that both have the same number of columns. On Thu, Oct 31, 2013 at 4:07 PM, Bert Gunter gunter.ber...@gene.com wrote: lapply() setdiff() by columns. Unless you have failed to tell us something, you almost certainly will not get a data frame (same number of rows/column) as your answer. -- Bert On Thu, Oct 31, 2013 at 12:58 PM, Yasin Gocgun entropy...@gmail.com wrote: Hi, I have two data frames, say, x and y, where y is a subset of x. How can I find the set difference of these two data frames (i.e., x-y)? Thanks, -- Yasin Gocgun __ 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. -- Bert Gunter Genentech Nonclinical Biostatistics (650) 467-7374 -- Yasin Gocgun __ 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-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] speed of makeCluster (package parallel)
See library(help = parallel”) On 28 Oct 2013, at 17:19, Arnaud Mosnier a.mosn...@gmail.com wrote: Hi all, I am quite new in the world of parallelization and I wonder if there is a way to increase the speed of creation of a parallel socket cluster. The time spend to include threads increase exponentially with the number of thread considered and I use of computer with two 8 cores CPU and thus showing a total of 32 threads in windows 7. Currently, I use the default parameters (type = PSOCK), but is there any fine tuning parameters that I can use to take advantage of this system ? Thanks in advance for your help ! Arnaud R version 3.0.1 (2013-05-16) Platform: x86_64-w64-mingw32/x64 (64-bit) [[alternative HTML version deleted]] __ 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-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] speed of makeCluster (package parallel)
First, use only the number of cores as a number of thread - i.e. I would not use hyper threading, etc.. Each core has its own caches and it is always fortunate if a process has enough memory; hyper threads use all the same cache on the core they are running on. detectCores() gives me for example 4 - I know I have 2. I would therefore call makeCluster() with nnode = 2. mcaffinity() lets you perform a technique called process-pinning (see process affinity) and is only possible if the OS supports it. It makes sometimes sense to assign certain processes to certain CPUs such that each process has enough memory in caches (e.g. for a 16 Core machine using 8 processes on CPUs 1, 3, 5, 7, 9, 11, 13 and 15; so each process has the cache of two CPUs). A lot of functions though are not available for Windows. At first it comes always the problem you want to solve and then you look how much memory will be used in a process and how much you have (more often the memory bandwidth is the bottleneck and not the computing power). Look at the architecture of your chips (how much L1 Cache, how much L2 cache). Then you decide how many cores to use and if it makes sense to pin processes to certain cores. There are no general recipes for parallel computing - each problem is different. Some problems are even not scalable. Simon On 28 Oct 2013, at 17:51, Arnaud Mosnier a.mosn...@gmail.com wrote: Thanks Simon, I already read the parallel vignette but I did not found what I wanted. May be you can be more specific on a part of the document that can provide me hints ! Arnaud 2013/10/28 Simon Zehnder szehn...@uni-bonn.de See library(help = parallel”) On 28 Oct 2013, at 17:19, Arnaud Mosnier a.mosn...@gmail.com wrote: Hi all, I am quite new in the world of parallelization and I wonder if there is a way to increase the speed of creation of a parallel socket cluster. The time spend to include threads increase exponentially with the number of thread considered and I use of computer with two 8 cores CPU and thus showing a total of 32 threads in windows 7. Currently, I use the default parameters (type = PSOCK), but is there any fine tuning parameters that I can use to take advantage of this system ? Thanks in advance for your help ! Arnaud R version 3.0.1 (2013-05-16) Platform: x86_64-w64-mingw32/x64 (64-bit) [[alternative HTML version deleted]] __ 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-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] DEoptim inconsistent output
It's always a good thing to trace the optimization. There can be a lot of reasons: Is the function correctly implemented? Is it defined for all values? If not, is it defined for all values in the constrained parameter space? Is it defined for the constraint? If it is not defined for the constrained (e.g. for zero, take the lower bound to be 1e-10; approximately zero) How do the numerical gradients behave? They rely on difference methods, could be undefined for some values. So, as long as we have no error messages or more precise information, we cannot help you. Best Simon On Sep 29, 2013, at 7:17 PM, Aya Anas aa...@feps.edu.eg wrote: Dear all, I wrote an R optimization code using the DEoptim package that is used to minimize a numerical integration subject to a single constraint. I noticed that the code takes a long time to give an output. The resulting output doesn't make sense at all. I got parameters that don't satisfy the constraint. In addition, when i substituted with the resulting parameters in the objective function, I didn't obtain the resulting optimal objective function value. What might be the reason behind this? . Thanks, Aya -- Best Regards Faculty of Economics Political Science Cairo University Tel:(202)35728055-(202)35728116-(202)35736608-(202)35736605 Fax:(202)35711020 Follow us on twitter:https://twitter.com/fepsnews [[alternative HTML version deleted]] __ 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-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] Memory distribution using foreach
First of all LSF is a batch scheduling software. It usually expects an .lsf script. Usually the compilers on a cluster are interchangeable via the 'module switch unload module load module' and MPI-2 is the message passing interface standard. This is also rather an topic for the high-performance R list. Next, doMC is a multicore package registering cores on one machine - AFAIK, i.e. you have to work on one machine with the 24 cores (inform yourself on the hardware on your cluster - there should also be introduction presentations online! To know what hardware you use and what architecture it has is the first step! Try 'bhosts' on your shell to see what hosts are available). If you want to use several machines, your backend for foreach should be doMPI and not doMC (see http://cran.r-project.org/web/packages/doMPI/vignettes/doMPI.pdf). If you found your host, you have to write an lsf-script like the one following (for OpenMP on ONE machine - using 24 cores, in most cases this suffices. Further, it is faster as you do not have to wait that long because you have to use just ONE machine. If you have BULL clusters - take these. A lot of cores 32/64… and a lot of memory) So in your case, write a script with extension .lsf containing: ### using the zsh shell #!/usr/bin/env zsh ### Job name #BSUB -J OpenMP ### File/path where output will be written, the %J is the job ID #BSUB -o OpenMP.%J ### (OFF) Different file for STDERR, if not to be merged with STDOUT # #BSUB -e OpenM.e%J ### Request the time you need for execution in minutes ### The format for the parameter is: [hour:]minute, ### that means for 80 minutes you could also use this: 1:20 #BSUB -W 3:00 ### Request virtual memory you need for your job in MB (per process) #BSUB -M 1024 ### Request higher amount of stack site (per process) #BSUB -S 1024 ### Request the number of compute slots you want to use #BSUB -n 24 ### Specify your mail address #BSUB pko...@bgc-jena.mpg.de ### Send a mail when job is done #BSUB -N ### Use esub for OpenMP #BSUB -a openmp ### (OFF) As R is usually compiled via gcc I would load the gcc module on your cluster # module switch pgi gcc/4.6 ### (OFF) load another OpenMP (check which one is usually loaded!! should be now OpenMP 4.0) version than the default one # module switch openmp openmp/3.0 ### Set stack and address limits ulimit -s unlimited ulimit -v unlimited ### Change to the work directory cd /home/your_username/ ### Execute your application (make sure, that R can be loaded via 'R' on the shell!!!) R --no-restore --no-save --quiet --slave your_R_script.R In your R script file, load the packages library(doMC) library(foreach) registerDoMC(24) ## now, foreach knows the backend. forach(...) %dopar% ….. ## save your stuff to your work- or home directory (csv or database) quit() --- Then you send the script to LSF via bsub - my_LSF_script.lsf Look via 'bjobs' if it is is send and what's its status (PEND or RUN). If the status is RUN you can look via 'bpeek your_job_ID' what the output looks like, while it runs. Best Simon On Sep 27, 2013, at 10:48 PM, pakoun pko...@bgc-jena.mpg.de wrote: Dear R users, I am struggling with memory issues and try to understand a few things. I am using an LSF cluster with PGI compiler and parallel mpi2 computing (whatever does that means..) and i submit a job like: bsub -R rusage[mem=3] -q queue -n 24 R CMD BATCH arguments.. myjob.r ..log According to that I am asking for 24 cores and 30GB RAM. Then I am using library(doMC) registerDoMC(24) and a foreach loop either simple or nested with the %dopar% command. 1. this 30 GB will be distributed among the 24 jobs or each will take 30? 2. If i dont ask the -n 24 argument still the foreach loop will run in parallel as i check with TOP command. What is the purpose of using it? Just to reserve the nodes from other users? Thank you -- View this message in context: http://r.789695.n4.nabble.com/Memory-distribution-using-foreach-tp4677133.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-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] Installing R, R packages
Had yesterday something pretty similar on the 73 installing benchmark: Error in untar2(tarfile, files, list, exdir, restore_times) : incomplete block on file The downloaded source packages are in ‘/private/var/folders/n9/zxfxcd01557dc06bf3c1njy0gn/T/RtmpkwQyuQ/downloaded_packages’ Warning messages: 1: In download.file(url, destfile, method, mode = wb, ...) : downloaded length 247786 != reported length 474551 2: In download.file(url, destfile, method, mode = wb, ...) : downloaded length 280942 != reported length 403808 3: In install.packages(benchmark) : installation of package ‘relations’ had non-zero exit status 4: In install.packages(benchmark) : installation of package ‘benchmark’ had non-zero exit status Today doing the same worked perfectly. On Sep 24, 2013, at 6:10 PM, David Winsemius dwinsem...@comcast.net wrote: On Sep 24, 2013, at 12:47 AM, Jeff Newmiller wrote: I have stopped using the Berkeley mirror, and just automatically use UCLA due to missing packages. However, I feel no compulsion to extrapolate and say that there is some sort of corruption going on at CRAN mirrors because it is only one data point. I have not been having this sort of difficulty with the Berkeley mirror. Never have seen this particular error. Occasionally see the mirror down on a weekend, but that happens with all mirrors. (I do use Berkeley as my standard mirror and I do download a relatively large number of packages especially at the point where I call `update.packages`.) -- David. Sent from my phone. Please excuse my brevity. David Arnold dwarnol...@suddenlink.net wrote: All, Consider this attempt: install.packages(car) trying URL 'http://cran.cnr.Berkeley.edu/bin/macosx/contrib/3.0/car_2.0-19.tgz' Content type 'application/x-gzip' length 1326903 bytes (1.3 Mb) opened URL == downloaded 1.1 Mb car/data/Rdata.rdb: Truncated tar archive tar: Error exit delayed from previous errors. The downloaded binary packages are in /var/folders/qE/qEavkZWTFMmxjncuY+HnqE+++TI/-Tmp-//Rtmpnxi1hV/downloaded_packages Warning messages: 1: In download.file(url, destfile, method, mode = wb, ...) : downloaded length 1126960 != reported length 1326903 2: 'tar' returned non-zero exit code 1 This seems to be happening frequently at the cran berkeley site. I've also had students try to install R from the Berkeley site and it just doesn't work for them. I had another student tell me today he tried all sorts of mirrors and could not get R and knitr installed until he tried the Washington site. Is there some sort of corruption going on at the CRAN mirrors? D. David Winsemius, MD Alameda, CA, USA __ 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-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] Checking a large data set for normality
Check the Jarque-Bera Test for univariate testing (http://hosho.ees.hokudai.ac.jp/~kubo/Rdoc/library/tseries/html/jarque.bera.test.html) and Mardia's test for multivariate testing (http://www.inside-r.org/packages/cran/MVN/docs/mardia.test). On Sep 24, 2013, at 6:44 PM, steric rd023...@ohio.edu wrote: Hello, I have a large data set that includes many soil parameters (i.e. pH, calcium levels, enzyme activity, etc) Does anyone have any input as to the easiest way to check a large data set for normality? Is there an R function/package that can do this all at once? Thank you in advance, -- View this message in context: http://r.789695.n4.nabble.com/Checking-a-large-data-set-for-normality-tp4676858.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-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] Hope u have some time for 2 more questions
Hi, as far as I know, there is no limitation on data size in regard to foreach. You should reserve though enough memory for your application on the cluster (via ulimit -s unlimited and ulimit -v unlimited). Furthermore I would check the following: Check if there are two versions of R on the cluster/your home directory on the frontend (LSF loads this frontend environment and uses the R version installed there). If you have two R executables (R and R64) make sure you use the 64bit version. Run R and call memory.limit() to see what are the limits of memory in your system. If this is limited to sizes below your needed sizes, increase it by calling R in the LSF script with the options --max-mem-size=YourSize and if you get errors of kind cannot allocate vector of size you should also use --max-vsize=YourVSize. Then, check if there is a memory leak in your application: If you compiled R with the --enable-memory-profiling you can use Rprof to do this otherwise you must rely on profiling instruments given by the cluster environment (I think you work there as well with modules, so type in the shell 'module avail' for listing available modules). If you detect a memory leak or if you see, that at certain points in your algorithms some objects are not used anymore call rm(ObjectName) and gc() for garbage collection. To your nested loop using foreach: That is a highly delicate issue in parallel computing and for the foreach syntax I refer to the must-read http://cran.r-project.org/web/packages/foreach/vignettes/nested.pdf. Using nested loops should be considered carefully in regard to organizing the nesting. In C++ you have the ability to determine how many cores should work on which loop. In the foreach environment using doMC this seems to me not possible. And, please keep the discussion to the r-help mailing list, so others can learn from it and researchers with more experience can also leave comments. Best Simon On Sep 19, 2013, at 9:24 PM, pko...@bgc-jena.mpg.de wrote: Hi again, if you have some time I would like to bother you again with 2 more questions. After your response the parallel code is working perfect but when I implement that to the real case (big matrices) I get an error for not numeric dimension and i guess that again it returns NULL or something. Are you aware if foreach loop can handle only a certain size objects? the equation that I am using includes 3 objects with 2Gb size each. The second question has to deal with the cores that foreach uses. Although I am asking to our cluster (LSF) to give me certain number of cpus, and also i am specifing that with library(doMC) registerDoMC(n) it seems from the top command that I am using all the cores. I am using 2 foreach as nest foreach(i in 1:16){ foreach(j in 1:10) etc etc.. maybe i should do something with this kind of nest? I am not aware about that. I am sorry for the long text , and thank you for your nice solution _ Sent from http://r.789695.n4.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] Renaming variables
You haven't said yet, what object your 'data file' is. If you mean a data.frame I would use colnames(dataName) - c(Col1Name, col2Name, ….) Best Simon On Sep 20, 2013, at 4:10 PM, Preetam Pal lordpree...@gmail.com wrote: Hi, I guess this is pretty basic. I have 25 variables in the data file (name: score), i.e. X1,X2,.,X25. I dont want to use score$X1, score$X2 everytime I use these variables. Is there a way I can rename all these variables as simply X1,X2,.X25 without writing 25 lines of code, one line for renaming each variable (eg: X1=score.X1 X2=score.X2 and so on) ? Thanks for your help. Regards, Preetam -- Preetam Pal (+91)-9432212774 M-Stat 2nd Year, Room No. N-114 Statistics Division, C.V.Raman Hall Indian Statistical Institute, B.H.O.S. Kolkata. [[alternative HTML version deleted]] __ 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-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] callNextMethod with dots argument
Kien, if you want to add variables in a function definition that is predefined by a Generic and calls CallNextMethod you have to add the '…' argument as well. setMethod(do,signature(a=Achild), function(a,msg,...) { print(do Achild) callNextMethod() }) Best Simon On Sep 19, 2013, at 8:58 AM, Kiên Kiêu kien.k...@jouy.inra.fr wrote: Hi, I met a problem when invoking callNextMethod within a method associated with a generic function taking ... as an argument. Here is the code setClass(Aparent,representation(x=numeric,y=numeric)) setClass(Achild,contains=Aparent) setGeneric(do,def=function(a,...) standardGeneric(do)) setMethod(do,signature(a=Aparent), function(a,msg) { print(do Aparent) }) setMethod(do,signature(a=Achild), function(a,msg) { print(do Achild) callNextMethod() }) myA - new(Achild) buf - do(a=myA) # works buf - do(a=myA,msg=bonjour) # error The last call yields the following error message: Error in callNextMethod() : in processing 'callNextMethod', found a '...' in the matched call, but no corresponding '...' argument which I do not understand. Replacing ... by msg in setGeneric makes it work. But I don't like this limitation so much (unless I understand it). Regards. Kien __ 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-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] R-3.0.1 g77 errors
On my systems Linux Scientific and Mac OS X I use as well for the F77 the gfortran compiler and this works. You could give it a trial. Best Simon On Sep 18, 2013, at 3:14 PM, Prigot, Jonathan jpri...@partners.org wrote: I am trying to build R-3.0.1 on our SPARC Solaris 10 system, but it fails part way through with g77 errors. Has anyone run into this? Any suggestions? For what it's worth, R-2.15.1 is the last one to build error free for us. === Jon Prigot R is now configured for sparc-sun-solaris2.10 Source directory: . Installation directory:/usr/local C compiler:gcc -std=gnu99 -g -O2 Fortran 77 compiler: g77 -g -O2 C++ compiler: g++ -g -O2 Fortran 90/95 compiler:gfortran Obj-C compiler: Interfaces supported: X11, tcltk External libraries:readline, ICU Additional capabilities: PNG, JPEG, TIFF, NLS Options enabled: shared BLAS, R profiling Recommended packages: yes make ... g77 -fPIC -g -O2 -ffloat-store -c dlamch.f -o dlamch.o dlamch.f: In function `dlamch': dlamch.f:89: warning: INTRINSIC DIGITS, EPSILON, HUGE, MAXEXPONENT, ^ Reference to unimplemented intrinsic `DIGITS' at (^) (assumed EXTERNAL) dlamch.f:89: INTRINSIC DIGITS, EPSILON, HUGE, MAXEXPONENT, ^ Invalid declaration of or reference to symbol `digits' at (^) [initially seen at (^)] dlamch.f:89: warning: INTRINSIC DIGITS, EPSILON, HUGE, MAXEXPONENT, ^ Reference to unimplemented intrinsic `EPSILON' at (^) (assumed EXTERNAL) -- Jonathan M. Prigot jpri...@partners.org Partners Healthcare Systems The information in this e-mail is intended only for th...{{dropped:18}} __ 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] foreach returns null first object in the list
Use an extra call to return: post.ls - foreach(i =1:2, .verbose = TRUE) %dopar% { fun - func.list[[i]] if (i == 1) return(fun(Xa, Sa)) if (i == 2) return(fun(Ta, Sa)) } Best Simon On Sep 17, 2013, at 11:55 AM, pakoun pko...@bgc-jena.mpg.de wrote: Dear all, I am sending you a copy of a demo that i am trying to make that work and embedded in a more complicated structure. The problem is that i want to get a list of 2 objects (matrices in specific, and thats why i am not specifing a .combine argument if i am not wrong..), but the first element of the list is just null. the second one is perfectly fine. Any idea? Thank you in advance library(doMC) registerDoMC(2) Xa-matrix(1:100,ncol=10) Sa-matrix(101:200,ncol=10) Ta-matrix(1:100,ncol=10) get.Xp - function(Xa,Sa){ result - Xa%*%Sa return(result) } get.Sa.post - function(Xa,Ta){ result - Ta%*%Sa return(result) } func.list - list(get.Xp,get.Sa.post) post.ls - foreach(i =1:2) %dopar% { fun - func.list[[i]] if(i==1){fun(Xa,Sa)} if(i==2){fun(Ta,Sa)} } post.ls [[1]] NULL [[2]] [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [1,] 154855 169455 184055 198655 213255 227855 242455 257055 271655 286255 [2,] 155910 170610 185310 200010 214710 229410 244110 258810 273510 288210 [3,] 156965 171765 186565 201365 216165 230965 245765 260565 275365 290165 [4,] 158020 172920 187820 202720 217620 232520 247420 262320 277220 292120 [5,] 159075 174075 189075 204075 219075 234075 249075 264075 279075 294075 [6,] 160130 175230 190330 205430 220530 235630 250730 265830 280930 296030 [7,] 161185 176385 191585 206785 221985 237185 252385 267585 282785 297985 [8,] 162240 177540 192840 208140 223440 238740 254040 269340 284640 299940 [9,] 163295 178695 194095 209495 224895 240295 255695 271095 286495 301895 [10,] 164350 179850 195350 210850 226350 241850 257350 272850 288350 303850 -- View this message in context: http://r.789695.n4.nabble.com/foreach-returns-null-first-object-in-the-list-tp4676303.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-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] foreach returns null first object in the list
Use an extra call to 'return()' as posted below. On Sep 17, 2013, at 5:55 PM, pakoun pko...@bgc-jena.mpg.de wrote: This is the output with the debugging message. I don't really understand what is the problem. post.ls - foreach(i =1:2, .verbose=T) %dopar% { + + fun - func.list[[i]] + if(i==1){fun(Xa,Sa)} + if(i==2){fun(Ta,Sa)} + + } numValues: 2, numResults: 0, stopped: TRUE got results for task 1 numValues: 2, numResults: 1, stopped: TRUE returning status FALSE got results for task 2 numValues: 2, numResults: 2, stopped: TRUE calling combine function evaluating call object to combine results: fun(accum, result.1, result.2) returning status TRUE -- View this message in context: http://r.789695.n4.nabble.com/foreach-returns-null-first-object-in-the-list-tp4676303p4676329.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-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] Unrecognized token
Maybe you should escape the single quotes something like ' id ' ? On Sep 17, 2013, at 6:03 PM, srecko joksimovic sreckojoksimo...@gmail.com wrote: Hi, when I generate query using sqldf library, like this: query = paste(paste(select * from tbl_user where student_id = , id, sep=), order by date_time, sep=) student - sqldf(query) everything works fine in case the id is 21328, 82882, or something like that. But, when id is something like 78789D, there is an error: Error in sqliteExecStatement(con, statement, bind.data) : RS-DBI driver: (error in statement: unrecognized token: 78789D) I tried replacing single quotes with double, but it still doesn't work... thanks, Srecko [[alternative HTML version deleted]] __ 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-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] set breakpoint in debug
You could just use debug(f) and then when the Browser opens and the loop begins type 'c', that jumps over the loop to next line after the loop. Best Simon On Sep 16, 2013, at 9:05 PM, Hui Du hui...@dataventures.com wrote: Hi All, I need some help regarding how to set up a breakpoint in debug. For example, I have a very simple/naïve function (a useless function just for demo) f = function() { x = 10; len = 100; a = 1; for(i in 1:len) { a = a * i; } y = x + a; y; } If I need to debug it, I can run debug(f). After I go into the debugger, if I want to skip the loop and stop in the statement y = y + a, directly, how to do that? I know R has a function named setBreakpoint but I have never used it correctly. Your help is highly appreciated. HXD [[alternative HTML version deleted]] __ 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-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] Package installation and path.package
I am following your suggestion and move this discussion to the R-devel list. Best Simon On Sep 9, 2013, at 7:58 AM, Prof Brian Ripley rip...@stats.ox.ac.uk wrote: On 09/09/2013 02:09, David Winsemius wrote: On Sep 8, 2013, at 8:00 AM, Simon Zehnder wrote: Dear R-Users and R-Devels, I am writing right now my own package that makes use of 'tempfile' and there within with 'path.package'. When I install it, I get the error: Error in path.package(mypackage) : none of the packages are loaded. I understand the error, but I would like to have a workaround. How can I give the path to the package I am actually installing without getting this error? (We do not have the code so this is speculation.) Your packages should be assumed to be available in one of the directories in .libPaths() Not until the package is fully installed. We do have no idea what is going on here ... and R-devel seems the appropriate list. -- Brian D. Ripley, rip...@stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ 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] Package installation and path.package
Dear R-Users and R-Devels, I am writing right now my own package that makes use of 'tempfile' and there within with 'path.package'. When I install it, I get the error: Error in path.package(mypackage) : none of the packages are loaded. I understand the error, but I would like to have a workaround. How can I give the path to the package I am actually installing without getting this error? Best Simon __ 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] [dfoptim] 'Error in fn(ginv(par), ...) : object 'alpha' not found'
Hi Carlos, your problem is a wrong definition of your Likelihood function. You call symbols in the code (alpha, beta) which have no value assigned to. When L the long calculation in the last lines is assigned to L alpha and beta do not exist. The code below corrects it. But you have a problem with a divergent integral when calling integrate. A problem you can surely fix as you know what your function is doing. Likelihood_cov - function(params, x, tx, T, IS) { r - params[1] alpha_zero - params[2] s - params[3] beta_zero - params[4] gamma_1 - params[5] gamma_2 - params[6] data$alpha - alpha_zero*exp(-gamma_1*IS) data$beta - beta_zero*exp(-gamma_2*IS) f - function(x, tx, T, alpha, beta) { g - function(y) (y + alpha)^(-( r + x))*(y + beta)^(-(s + 1)) integrate(g, tx, T)$value } integral - mdply(data, f) L - exp(lgamma(r+x)-lgamma(r)+r*(log(alpha_zero)-log(alpha_zero+T))-x*log(alpha_zero+T)+s*(log(beta_zero)-log(beta_zero+T)))+exp(lgamma(r+x)-lgamma(r)+r*log(alpha_zero)+log(s)+s*log(beta_zero)+log(integral$V1)) f - -sum(log(L)) return (f) } Best Simon On Sep 3, 2013, at 1:28 PM, Carlos Nasher carlos.nas...@googlemail.com wrote: Dear R helpers, I have problems to properly define a Likelihood function. Thanks to your help my basic model is running quite well, but I have problems to get the enhanced version (now incorporating covariates) running. Within my likelihood function I define a variable 'alpha'. When I want to optimize the function I get the error message: 'Error in fn(ginv(par), ...) : object 'alpha' not found' I think it's actually not a problem with the optimization function (nmkb), but with the Likelihood function itself. I do not understand why 'alpha' is a missing object. 'alpha' should be part of the dataframe 'data' (as 'beta' should be too), like 'x', 'tx', ''T. But it obviously isn't. Here's a minimum example which reproduces my problem: ## library(plyr) library(dfoptim) ### Sample data ### x - c(3, 0, 2, 5, 1, 0, 0, 1, 0, 2) tx - c(24.57, 0.00, 26.86, 34.57, 2.14, 0.00, 0.00, 8.57, 0.00, 14.29) T - c(33.29, 30.71, 31.29, 34.57, 36.00, 35.43, 31.14, 33.86, 35.71, 35.86) IS - c(54.97, 13.97, 122.33, 110.84, 30.72, 14.96, 30.72, 20.74, 29.16, 83.00) data - data.frame(x=x, tx=tx, T=T) rm(x, tx, T) ### Likelihood function ### Likelihood_cov - function(params, x, tx, T, IS) { r - params[1] alpha_zero - params[2] s - params[3] beta_zero - params[4] gamma_1 - params[5] gamma_2 - params[6] data$alpha - alpha_zero*exp(-gamma_1*IS) data$beta - beta_zero*exp(-gamma_2*IS) f - function(x, tx, T, alpha, beta) { g - function(y) (y + alpha)^(-( r + x))*(y + beta)^(-(s + 1)) integrate(g, tx, T)$value } integral - mdply(data, f) L - exp(lgamma(r+x)-lgamma(r)+r*(log(alpha)-log(alpha+T))-x*log(alpha+T)+s*(log(beta)-log(beta+T)))+exp(lgamma(r+x)-lgamma(r)+r*log(alpha)+log(s)+s*log(beta)+log(integral$V1)) f - -sum(log(L)) return (f) } ### ML optimization ### params - c(0.2, 5, 0.2, 5, -0.02, -0.02) fit - nmkb(par=params, fn=Likelihood_cov, lower=c(0.0001, 0.0001, 0.0001, 0.0001, -Inf, -Inf), upper=c(Inf, Inf, Inf, Inf, Inf, Inf), x=data$x, tx=data$tx, T=data$T, IS=IS) ## Maybe you could give me a hint were the flaw in my code is. Many thanks in advance. Carlos - Carlos Nasher Buchenstr. 12 22299 Hamburg tel:+49 (0)40 67952962 mobil:+49 (0)175 9386725 mail: carlos.nas...@gmail.com [[alternative HTML version deleted]] __ 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-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] How to catch errors regarding the hessian in 'optim'
Dear John, thank you very much for your answer. I take a look at these packages (Rvmmin and Rcgmin). That sounds very interesting. For the example: The method relies on data which I always try to avoid to send on the r-help list - not that my data is confidential - but it becomes even more cumbersome when sending datasets over the list. I made the experience, that in such case answers are rare. Maybe you have a suggestion if and how users should send data with their code. I am interested in your opinion and thankful for sharing it. Best Simon On Sep 2, 2013, at 4:42 PM, Prof J C Nash (U30A) nas...@uottawa.ca wrote: This may be one of the many mysteries of the internals of L-BFGS-B, which I have found fails from time to time. That is one of the reasons for Rvmmin and Rcgmin (and hopefully sooner rather than later Rtn - a truncated Newton method, currently working for unconstrained problems, but still glitchy for bounds constraints). These are all-R codes so that users and developers can get inside to control special situations. If you have a test problem -- the infamous reproducible example -- there are several of us who can likely help to sort out your troubles. JN On 13-09-02 06:00 AM, r-help-requ...@r-project.org wrote: Message: 10 Date: Sun, 1 Sep 2013 17:09:35 +0200 From: Simon Zehnderszehn...@uni-bonn.de To: R-help helpr-help@r-project.org Subject: [R] How to catch errors regarding the hessian in 'optim' Message-ID:eb37670e-8544-4c89-9172-245eb6cc5...@uni-bonn.de Content-Type: text/plain; charset=us-ascii Dear R-Users and R-Developers, in a comparison between two different estimation approaches I would like to catch errors from optim regarding the hessian matrix. I use optim with method = L-BFGS-B thereby relying on numerical differentiation for the hessian matrix. I do know, that the estimation approach that uses numerical optimization has sometimes problems with singular hessian matrices and I consider it as one of its disadvantages of this method. To show the frequency of such problems in my simulation study I have to set 'hessian = TRUE' and to collect the errors from optim regarding the hessian. Now I am a little stucked how I could catch specifically errors from the hessian matrix in 'optim'. I do know that such errors are thrown most certainly from function 'La_solve' in Lapack.c. Does anyone has an idea how I could solve this task (clearly with tryCatch but how to choose only errors for the hessian)? Best Simon __ 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] How to catch errors regarding the hessian in 'optim'
Dear R-Users and R-Developers, in a comparison between two different estimation approaches I would like to catch errors from optim regarding the hessian matrix. I use optim with method = L-BFGS-B thereby relying on numerical differentiation for the hessian matrix. I do know, that the estimation approach that uses numerical optimization has sometimes problems with singular hessian matrices and I consider it as one of its disadvantages of this method. To show the frequency of such problems in my simulation study I have to set 'hessian = TRUE' and to collect the errors from optim regarding the hessian. Now I am a little stucked how I could catch specifically errors from the hessian matrix in 'optim'. I do know that such errors are thrown most certainly from function 'La_solve' in Lapack.c. Does anyone has an idea how I could solve this task (clearly with tryCatch but how to choose only errors for the hessian)? Best Simon __ 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] Understanding S4 method dispatch
Ambiguity is indeed detected by R and the user is informed on it. But in the case of Hadley's example, I still believe, that the specific multiple inheritance structure creates this behavior. If you call: showMethods(f) Function: f (package .GlobalEnv) x=A, y=A x=AB, y=AB (inherited from: x=B, y=B) x=B, y=B you see, that AB inherited the method from B. As B inherits already from A and AB does as well, AB would have two A objects, the one from its inheritance from A and the one from B. In the case of direct inheritance, the only slot is the .xData slot, which is NULL. If one class contains another one, it contains all its slots. As the .xData slot is the same for each class A, B and AB, R must remove the ambiguity of slots, this is done by allowing only one slot .xData, which is the one from B. So in some way A is hidden by B in AB - but this is actually a common technique in OOP. In C++ for example we have the same problem and we can define which members should be inherited if multiple inheritance is aware, by using the keyword 'virtual'. Towards Hervé's second point: If we build a class C: setClass(C, contains = c(A)) setMethod(f, signature(C, C), function(x, y) C-C) And then construct a multiple inheritance structure within a class BC: setClass(BC, contains = c(B, C)) bc - new(BC) we see, that indeed the first lexicographical signature is chosen: showMethods(f) Function: f (package .GlobalEnv) x=A, y=A x=B, y=B x=BC, y=BC (inherited from: x=B, y=B) x=C, y=C f(bc, bc) [1] B-B In this case, we have a 'true' tie in the inheritance structure: B from A and C from A, any one of the two As have to be taken into the BC object and the one from B is chosen. But method dispatch is in this case independent of the second level inheritance. B and C are not related directly to each other, so method dispatch is chosen lexicographically. In my opinion the reason for the behavior lies in the specific multiple inheritance structure between AB, B and A. Best Simon On Aug 14, 2013, at 2:11 AM, Hervé Pagès hpa...@fhcrc.org wrote: Hi Hadley, I suspect that the dispatch algorithm doesn't realize that selection is ambiguous in your example. For 2 reasons: (1) When it does realize it, it notifies the user: setClass(A, NULL) setGeneric(f, function(x, y) standardGeneric(f)) setMethod(f, signature(A, ANY), function(x, y) A-ANY) setMethod(f, signature(ANY, A), function(x, y) ANY-A) a - new(A) Then: f(a, a) Note: method with signature ‘A#ANY’ chosen for function ‘f’, target signature ‘A#A’. ANY#A would also be valid [1] A-ANY (2) When dispatch is ambiguous, the first method lexicographically in the ordering should be selected (according to ?Methods). So it should be A#A, not B#B. So it looks like a bug to me... Cheers, H. On 08/13/2013 06:08 AM, Hadley Wickham wrote: Hi all, Any insight into the code below would be appreciated - I don't understand why two methods which I think should have equal distance from the call don't. Thanks! Hadley # Create simple class hierarchy setClass(A, NULL) setClass(B, A) a - new(A) b - new(B) setGeneric(f, function(x, y) standardGeneric(f)) setMethod(f, signature(A, A), function(x, y) A-A) setMethod(f, signature(B, B), function(x, y) B-B) # These work as I expect f(a, a) f(b, b) setClass(AB, contains = c(A, B)) ab - new(AB) # Why does this return B-B? Shouldn't both methods be an equal distance? f(ab, ab) # These both return distance 1, as I expected extends(AB, A, fullInfo=TRUE)@distance extends(AB, B, fullInfo=TRUE)@distance # So why is signature(B, B) closer than signature(A, A) -- Hervé Pagès Program in Computational Biology Division of Public Health Sciences Fred Hutchinson Cancer Research Center 1100 Fairview Ave. N, M1-B514 P.O. Box 19024 Seattle, WA 98109-1024 E-mail: hpa...@fhcrc.org Phone: (206) 667-5791 Fax:(206) 667-1319 __ 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-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] Grabbing numbers inside a character string
Thomas, se ?sub and regular expression. That should make it. Further, see the package gsubfn Best Simon On Aug 14, 2013, at 2:14 PM, Thomas thomas.ches...@nottingham.ac.uk wrote: I have a string that contains something like: ...verified email at neu.edubrCited by 99853br/td/tr/table/div/div/div/divdiv... and I'd like to extract the number next to the text Cited by - so it will be whatever numbers are beside Cited by until a non-numeric character is reached. Can anyone suggest a good way of doing this please? Thank you, Thomas This message and any attachment are intended solely for the addressee and may contain confidential information. If you have received this message in error, please send it back to me, and immediately delete it. Please do not use, copy or disclose the information contained in this message or in any attachment. Any views or opinions expressed by the author of this email do not necessarily reflect the views of the University of Nottingham. This message has been checked for viruses but the contents of an attachment may still contain software viruses which could damage your computer system, you are advised to perform your own checks. Email communications with the University of Nottingham may be monitored as permitted by UK legislation. __ 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-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] Understanding S4 method dispatch
Because the signature is always (A,A) or (B,B). Then, as in AB we have A and B and no relationship between A and B, R chooses the method lexicographically. The result is as expected: f for A is chosen. If you would do something like: setClass(A, contains = list) setClass(B, contains = list) setClass(AB, contains = c(A, B)) setGeneric(f, function(x, y) standardGeneric(f)) setMethod(f, signature(A, ANY), function(x, y) A-ANY) setMethod(f, signature(ANY, A), function(x, y) ANY-A) setMethod(f, signature(B, B), function(x, y) B-B) ab - new(AB) f(ab, ab) You get an ambiguity, as there is no function with signature pair that can be directly matched for A. But for B such a signature pair exists. So R chooses B. If you would specify again a function for signature (A,A) we are back: setMethod(f, signature(A, A), function(x, y) A-A) f(ab, ab) Best Simon On Aug 14, 2013, at 5:25 PM, Hadley Wickham h.wick...@gmail.com wrote: In my opinion the reason for the behavior lies in the specific multiple inheritance structure between AB, B and A. So what if we don't make such a weird inheritance structure, and instead have A and B inherit from a common parent: setClass(A, contains = list) setClass(B, contains = list) setClass(AB, contains = c(A, B)) setGeneric(f, function(x, y) standardGeneric(f)) setMethod(f, signature(A, A), function(x, y) A-A) setMethod(f, signature(B, B), function(x, y) B-B) ab - new(AB) f(ab, ab) Why isn't there a warning about ambiguous dispatch? Hadley -- Chief Scientist, RStudio http://had.co.nz/ __ 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] Understanding S4 method dispatch
Hadley, The class AB inherits from A and from B, but B already inherits from class A. So actually you only have an object of class B in your object of class AB. When you call the function f R looks for a method f for AB objects. It does not find such a method and looks for a method of the object inherited from, B. Such a method is present and is then executed. The inheritance structure has to be changed. The behavior is actually desired, as if this behavior weren't given a diamond class inheritance would be fatal. Best Simon On Aug 13, 2013, at 3:08 PM, Hadley Wickham h.wick...@gmail.com wrote: Hi all, Any insight into the code below would be appreciated - I don't understand why two methods which I think should have equal distance from the call don't. Thanks! Hadley # Create simple class hierarchy setClass(A, NULL) setClass(B, A) a - new(A) b - new(B) setGeneric(f, function(x, y) standardGeneric(f)) setMethod(f, signature(A, A), function(x, y) A-A) setMethod(f, signature(B, B), function(x, y) B-B) # These work as I expect f(a, a) f(b, b) setClass(AB, contains = c(A, B)) ab - new(AB) # Why does this return B-B? Shouldn't both methods be an equal distance? f(ab, ab) # These both return distance 1, as I expected extends(AB, A, fullInfo=TRUE)@distance extends(AB, B, fullInfo=TRUE)@distance # So why is signature(B, B) closer than signature(A, A) -- Chief Scientist, RStudio http://had.co.nz/ __ 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-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] Understanding S4 method dispatch
If you take an example which works with slots, setClass(A, representation(a = numeric) setClass(B, contains = c(A), representation(b = numeric)) a - new(A, a = 2) b - new(B, a = 3, b = 2) setClass(AB, contains = c(A, B)) new(AB, a = 2, b = 3) You see, that there is only one @a slot, the one inherited from B, that B inherits from A. If this were not the case, which slot should be taken, if we would call @a? To avoid this kind of ambiguity, only one A class is inherited to AB: the one B already inherits from A. You could create a class, that contains another A object in a slot: setClass(AandB, contains = c(B), representation(A = A)) new(AandB, a = 2, b = 3, A = new(A, a = 3)) Now back to your example: as there is only one A object inside the B object which is contained by the AB object, method dispatch works the way as it should: It looks for a method f for an AB object. It does not find one. Then it looks for a method f for the contained B object (as this is the only one contained in AB) and it finds a method. Then it calls this method on the B part of the object AB and the result is B-B Best Simon On Aug 13, 2013, at 4:24 PM, Hadley Wickham h.wick...@gmail.com wrote: The class AB inherits from A and from B, but B already inherits from class A. So actually you only have an object of class B in your object of class AB. When you call the function f R looks for a method f for AB objects. It does not find such a method and looks for a method of the object inherited from, B. Such a method is present and is then executed. The inheritance structure has to be changed. The behavior is actually desired, as if this behavior weren't given a diamond class inheritance would be fatal. Are you sure? That behaviour doesn't agree with the description of method dispatch given in ?Methods, not with getClass(AB) which shows that AB inherits from both A and B. (I totally agree that this is a bad idea, and unlikely to be useful in real life, but I'm trying to understand the details of S4 dispatch) getClass(AB) Class AB [in .GlobalEnv] Slots: Name: .xData Class: NULL Extends: Class B, directly Class A, directly Class .NULL, by class A, distance 2 Class NULL, by class A, distance 3, with explicit coerce Class OptionalFunction, by class A, distance 4, with explicit coerce Class optionalMethod, by class A, distance 4, with explicit coerce -- Chief Scientist, RStudio http://had.co.nz/ __ 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] Outliers and overdispersion
I do not know what you are exactly estimating, but if it is about count models and the model fit gets better when you drop the outliers, it does not say, that the model is now more correct. It just says, if the data were without the outliers, this model would fit good. Overdispersion in count data is sometimes a cue, that you have a mixture distribution as the generating process - for example instead of one, K different (sub)species of birds which were aggregated in the count data. In this case a mixture (negative binomial)- distribution with K components could fit the data better. Best Simon On Aug 13, 2013, at 5:28 PM, Marta Lomas lomasv...@hotmail.com wrote: Hi again, I have a question on some outliers that I have in my response variable (wich are bird counts). At the beginning I did not drop them out because they are part of the normal counts and I considered them ecologically correct. However, I tried some of the same models without ouliers and the AICs are thus better. I also have nice significances this way... So would you say that, even though the outliers are right observations and taking into consideration that already the negative binomial distribution that I am using is accounting for the some of the overdispersion due to the outliers, it is better to drop them out as the models fit better this way? Thanks for your patience! :) [[alternative HTML version deleted]] __ 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-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] decimal separator from comma to dot
Hi, I think this could help you: https://stat.ethz.ch/pipermail/r-help/2008-January/152378.html Best Simon On Aug 9, 2013, at 12:19 PM, maxbre mbres...@arpa.veneto.it wrote: This is my reproducible example df-structure(list(IDANT = c(37837L, 37838L, 37839L, 37840L, 37841L, 37842L, 37843L, 40720L, 40721L, 40722L), N_TX = c(6L, 6L, 6L, 4L, 1L, 1L, 1L, 2L, 2L, 1L), TILT = c(0L, 0L, 0L, 0L, 6L, 6L, 6L, 0L, 0L, 0L), DIREZIONE = c(50L, 220L, 110L, 50L, 220L, 110L, 50L, 170L, 70L, 270L), DATA_INI = structure(c(2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L), .Label = c(20/10/2004, 29/08/2002 ), class = factor), POT_TX = structure(c(4L, 4L, 4L, 3L, 2L, 2L, 2L, 1L, 1L, 1L), .Label = c(10, 11,5, 4, 8), class = factor)), .Names = c(IDANT, N_TX, TILT, DIREZIONE, DATA_INI, POT_TX), row.names = c(NA, 10L), class = data.frame) The data frame “df” it’s a “simplified snapshot” from an oracle table imported via odbc (by using the package RODBC) ; now, the problem with the data frame “df” is that one variable “POT_TX” have the decimal separator formatted with comma instead of dot; so I’m in the need to reformat the variable from factor to numeric type to perform some useful calculations this is my code I worked so far #function to change comma to dot myfun - function(x) {sub(,,.,x)} #apply the function to all variables new - apply(df, 2, myfun ) newdf - data.frame(apply(new, 2, as.numeric)) str(newdf) #apply the function to one variables var-as.numeric(sapply(df[,6], FUN=myfun)) df$POT_TX-var str(df) my questions: 1.is it possible to use the more convenient function “apply” with reference to just some variables and not all of them? I’ve been reading on the help that “Where X has named dimnames, it can be a character vector selecting dimension names.”, what does exactly means that? 2.do you know some more convenient methods to perform comma to dot substitution in just some selected variables of a data frame? To note that the desired final result is again a dataframe 3.do you know any workaround to set the decimal delimiter for a variable (field) “at the early steps” when connecting via odbc and then querying to the table? For these “early steps” I’ve been using something like: library(RODBC) con-odbcConnect(dsn, uid=user, pwd=password) df-sqlQuery(con, select * from table.name;) thank you all for the eventual support -- View this message in context: http://r.789695.n4.nabble.com/decimal-separator-from-comma-to-dot-tp4673414.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-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] Method dispatch in S4
Hi Martin, thank you very much for this profound answer! Your added design advice is very helpful, too! For the 'simple example': Sometimes I am still a little overwhelmed from a certain setting in the code and my ideas how I want to handle a process. But I learn from session to session. In future I will also span the lines more than 80 columns. I am used to the indent in my vim editor. I have one further issue: I do know, that you are one of the leading developers of the bioconductor package which uses (as far as I have read) extensively OOP in R. Is there a package you could suggest to me to learn from by reading and understanding the code? Where can I find the source code? Best Simon On Aug 8, 2013, at 10:00 PM, Martin Morgan mtmor...@fhcrc.org wrote: On 08/04/2013 02:13 AM, Simon Zehnder wrote: So, I found a solution: First in the initialize method of class C coerce the C object into a B object. Then call the next method in the list with the B class object. Now, in the initialize method of class B the object is a B object and the respective generateSpec method is called. Then, in the initialize method of C the returned object from callNextMethod has to be written to the C class object in .Object. See the code below. setMethod(initialize, C, function(.Object, value) {.Object@c - value; object - as(.Object, B); object - callNextMethod(object, value); as(.Object, B) - object; .Object - generateSpec(.Object); return(.Object)}) This setting works. I do not know though, if this setting is the usual way such things are done in R OOP. Maybe the whole class design is disadvantageous. If anyone detects a mistaken design, I am very thankful to learn. Hi Simon -- your 'simple' example is pretty complicated, and I didn't really follow it in detail! The code is not formatted for easy reading (e.g., lines spanning no more than 80 columns) and some of it (e.g., generateSpec) might not be necessary to describe the problem you're having. A good strategy is to ensure that 'new' called with no arguments works (there are other solutions, but following this rule has helped me to keep my classes and methods simple). This is not the case for new(A) new(C) The reason for this strategy has to do with the way inheritance is implemented, in particular the coercion from derived to super class. Usually it is better to provide default values for arguments to initialize, and to specify arguments after a '...'. This means that your initialize methods will respects the contract set out in ?initialize, in particular the handling of unnamed arguments: ...: data to include in the new object. Named arguments correspond to slots in the class definition. Unnamed arguments must be objects from classes that this class extends. I might have written initialize,A-method as setMethod(initialize, A, function(.Object, ..., value=numeric()){ .Object - callNextMethod(.Object, ..., a=value) generateSpec(.Object) }) Likely in a subsequent iteration I would have ended up with (using the convention that function names preceded by '.' are not exported) .A - setClass(A, representation(a = numeric, specA = numeric)) .generateSpecA - function(a) { 1 / a } A - function(a=numeric(), ...) { specA - .generateSpecA(a) .A(..., a=a, specA=specA) } setMethod(generateSpec, A, function(object) { .generateSpecA(object@a) }) ensuring that A() returns a valid object and avoiding the definition of an initialize method entirely. Martin Best Simon On Aug 3, 2013, at 9:43 PM, Simon Zehnder simon.zehn...@googlemail.com wrote: setMethod(initialize, C, function(.Object, value) {.Object@c - value; .Object - callNextMethod(.Object, value); .Object - generateSpec(.Object); return(.Object)}) __ 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. -- Computational Biology / Fred Hutchinson Cancer Research Center 1100 Fairview Ave. N. PO Box 19024 Seattle, WA 98109 Location: Arnold Building M1 B861 Phone: (206) 667-2793 __ 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] Method dispatch in S4
Hi Bert, thank you very much for your suggestion! I will take a look at it soon! Best Simon On Aug 9, 2013, at 4:45 PM, Bert Gunter gunter.ber...@gene.com wrote: Simon: Have a look at the proto package for which there is a vignette. You may find it suitable for your needs and less intimidating. Cheers, Bert On Fri, Aug 9, 2013 at 7:40 AM, Simon Zehnder szehn...@uni-bonn.de wrote: Hi Martin, thank you very much for this profound answer! Your added design advice is very helpful, too! For the 'simple example': Sometimes I am still a little overwhelmed from a certain setting in the code and my ideas how I want to handle a process. But I learn from session to session. In future I will also span the lines more than 80 columns. I am used to the indent in my vim editor. I have one further issue: I do know, that you are one of the leading developers of the bioconductor package which uses (as far as I have read) extensively OOP in R. Is there a package you could suggest to me to learn from by reading and understanding the code? Where can I find the source code? Best Simon On Aug 8, 2013, at 10:00 PM, Martin Morgan mtmor...@fhcrc.org wrote: On 08/04/2013 02:13 AM, Simon Zehnder wrote: So, I found a solution: First in the initialize method of class C coerce the C object into a B object. Then call the next method in the list with the B class object. Now, in the initialize method of class B the object is a B object and the respective generateSpec method is called. Then, in the initialize method of C the returned object from callNextMethod has to be written to the C class object in .Object. See the code below. setMethod(initialize, C, function(.Object, value) {.Object@c - value; object - as(.Object, B); object - callNextMethod(object, value); as(.Object, B) - object; .Object - generateSpec(.Object); return(.Object)}) This setting works. I do not know though, if this setting is the usual way such things are done in R OOP. Maybe the whole class design is disadvantageous. If anyone detects a mistaken design, I am very thankful to learn. Hi Simon -- your 'simple' example is pretty complicated, and I didn't really follow it in detail! The code is not formatted for easy reading (e.g., lines spanning no more than 80 columns) and some of it (e.g., generateSpec) might not be necessary to describe the problem you're having. A good strategy is to ensure that 'new' called with no arguments works (there are other solutions, but following this rule has helped me to keep my classes and methods simple). This is not the case for new(A) new(C) The reason for this strategy has to do with the way inheritance is implemented, in particular the coercion from derived to super class. Usually it is better to provide default values for arguments to initialize, and to specify arguments after a '...'. This means that your initialize methods will respects the contract set out in ?initialize, in particular the handling of unnamed arguments: ...: data to include in the new object. Named arguments correspond to slots in the class definition. Unnamed arguments must be objects from classes that this class extends. I might have written initialize,A-method as setMethod(initialize, A, function(.Object, ..., value=numeric()){ .Object - callNextMethod(.Object, ..., a=value) generateSpec(.Object) }) Likely in a subsequent iteration I would have ended up with (using the convention that function names preceded by '.' are not exported) .A - setClass(A, representation(a = numeric, specA = numeric)) .generateSpecA - function(a) { 1 / a } A - function(a=numeric(), ...) { specA - .generateSpecA(a) .A(..., a=a, specA=specA) } setMethod(generateSpec, A, function(object) { .generateSpecA(object@a) }) ensuring that A() returns a valid object and avoiding the definition of an initialize method entirely. Martin Best Simon On Aug 3, 2013, at 9:43 PM, Simon Zehnder simon.zehn...@googlemail.com wrote: setMethod(initialize, C, function(.Object, value) {.Object@c - value; .Object - callNextMethod(.Object, value); .Object - generateSpec(.Object); return(.Object)}) __ 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. -- Computational Biology / Fred Hutchinson Cancer Research Center 1100 Fairview Ave. N. PO Box 19024 Seattle, WA 98109 Location: Arnold Building M1 B861 Phone: (206) 667-2793 __ 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
Re: [R] Method dispatch in S4
Hi Martin, is proto in S3? I will take a look first at the simple package EBImage. Thank you very much for the suggestions! Best Simon On Aug 9, 2013, at 5:01 PM, Martin Morgan mtmor...@fhcrc.org wrote: On 08/09/2013 07:45 AM, Bert Gunter wrote: Simon: Have a look at the proto package for which there is a vignette. You may find it suitable for your needs and less intimidating. Won't help much with S4, though! Some answers here http://stackoverflow.com/questions/5437238/which-packages-make-good-use-of-s4-objects including from Bioconductor simple class in EBImage, the advanced IRanges package and the 'toy' StudentGWAS. Martin Cheers, Bert On Fri, Aug 9, 2013 at 7:40 AM, Simon Zehnder szehn...@uni-bonn.de wrote: Hi Martin, thank you very much for this profound answer! Your added design advice is very helpful, too! For the 'simple example': Sometimes I am still a little overwhelmed from a certain setting in the code and my ideas how I want to handle a process. But I learn from session to session. In future I will also span the lines more than 80 columns. I am used to the indent in my vim editor. I have one further issue: I do know, that you are one of the leading developers of the bioconductor package which uses (as far as I have read) extensively OOP in R. Is there a package you could suggest to me to learn from by reading and understanding the code? Where can I find the source code? Best Simon On Aug 8, 2013, at 10:00 PM, Martin Morgan mtmor...@fhcrc.org wrote: On 08/04/2013 02:13 AM, Simon Zehnder wrote: So, I found a solution: First in the initialize method of class C coerce the C object into a B object. Then call the next method in the list with the B class object. Now, in the initialize method of class B the object is a B object and the respective generateSpec method is called. Then, in the initialize method of C the returned object from callNextMethod has to be written to the C class object in .Object. See the code below. setMethod(initialize, C, function(.Object, value) {.Object@c - value; object - as(.Object, B); object - callNextMethod(object, value); as(.Object, B) - object; .Object - generateSpec(.Object); return(.Object)}) This setting works. I do not know though, if this setting is the usual way such things are done in R OOP. Maybe the whole class design is disadvantageous. If anyone detects a mistaken design, I am very thankful to learn. Hi Simon -- your 'simple' example is pretty complicated, and I didn't really follow it in detail! The code is not formatted for easy reading (e.g., lines spanning no more than 80 columns) and some of it (e.g., generateSpec) might not be necessary to describe the problem you're having. A good strategy is to ensure that 'new' called with no arguments works (there are other solutions, but following this rule has helped me to keep my classes and methods simple). This is not the case for new(A) new(C) The reason for this strategy has to do with the way inheritance is implemented, in particular the coercion from derived to super class. Usually it is better to provide default values for arguments to initialize, and to specify arguments after a '...'. This means that your initialize methods will respects the contract set out in ?initialize, in particular the handling of unnamed arguments: ...: data to include in the new object. Named arguments correspond to slots in the class definition. Unnamed arguments must be objects from classes that this class extends. I might have written initialize,A-method as setMethod(initialize, A, function(.Object, ..., value=numeric()){ .Object - callNextMethod(.Object, ..., a=value) generateSpec(.Object) }) Likely in a subsequent iteration I would have ended up with (using the convention that function names preceded by '.' are not exported) .A - setClass(A, representation(a = numeric, specA = numeric)) .generateSpecA - function(a) { 1 / a } A - function(a=numeric(), ...) { specA - .generateSpecA(a) .A(..., a=a, specA=specA) } setMethod(generateSpec, A, function(object) { .generateSpecA(object@a) }) ensuring that A() returns a valid object and avoiding the definition of an initialize method entirely. Martin Best Simon On Aug 3, 2013, at 9:43 PM, Simon Zehnder simon.zehn...@googlemail.com wrote: setMethod(initialize, C, function(.Object, value) {.Object@c - value; .Object - callNextMethod(.Object, value); .Object - generateSpec(.Object); return(.Object)}) __ 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. -- Computational Biology
Re: [R] Method dispatch in S4
So, I found a solution: First in the initialize method of class C coerce the C object into a B object. Then call the next method in the list with the B class object. Now, in the initialize method of class B the object is a B object and the respective generateSpec method is called. Then, in the initialize method of C the returned object from callNextMethod has to be written to the C class object in .Object. See the code below. setMethod(initialize, C, function(.Object, value) {.Object@c - value; object - as(.Object, B); object - callNextMethod(object, value); as(.Object, B) - object; .Object - generateSpec(.Object); return(.Object)}) This setting works. I do not know though, if this setting is the usual way such things are done in R OOP. Maybe the whole class design is disadvantageous. If anyone detects a mistaken design, I am very thankful to learn. Best Simon On Aug 3, 2013, at 9:43 PM, Simon Zehnder simon.zehn...@googlemail.com wrote: setMethod(initialize, C, function(.Object, value) {.Object@c - value; .Object - callNextMethod(.Object, value); .Object - generateSpec(.Object); return(.Object)}) __ 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] fortran random number
Hi Filippo, I would advice you to use the RcppArmadillo package - maybe in combination with the inline package. It is extremely fast and very flexible. A simple RNGScope() lets you take the RNG from R. I myself programmed an MCMC sampler with this package and it is much faster, than the one I had before in matlab. Rcpp gives you a very nice API to communicate with R + you can reuse memory from R which makes it even faster. You find a lot of documentation online: http://www.rcpp.org http://cran.r-project.org/web/packages/RcppArmadillo/index.html http://cran.r-project.org/web/packages/Rcpp/index.html Best Simon On Aug 4, 2013, at 10:38 PM, Filippo ingf...@gmail.com wrote: Hi, I need an advice. I'm writing a Markov Chain Monte Carlo subroutine in F95 which should be called from R. I noticed that the subroutine RANDOM_NUMBER takes the seed from R (is that correct?). I was wondering if is better initialize the random seed inside the F95 subroutine or use the R seed. Regards, Filippo __ 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-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] Method dispatch in S4
Dear R-Users and R-Devels, I am struggling with the method dispatch for S4 objects. First I will give a simple example to make clear what the general setting in my project is. Three classes will be build in the example: A is just a simple class. B contains A in a slot and C inherits from B. In addition all classes have a slot called specA, specB, specC that should be computed during initialization by calling the method generateSpec defined for each class A,B,C: ## Define class A## setClass(A, representation(a = numeric, specA = numeric)) setMethod(initialize, A, function(.Object, value){.Object@a - value; .Object - generateSpec(.Object); return(.Object)}) setGeneric(generateSpec, function(object) standardGeneric(generateSpec)) setMethod(generateSpec, A, function(object) {object@specA - 1/object@a; return(object)}) ## Define class B containing A in a slot ## setClass(B, representation(b = numeric, specB = numeric, A = A)) setMethod(initialize, B, function(.Object, value = 0){.Object@b - value; .Object@A - new(A, value = value); .Object - generateSpec(.Object); return(.Object)}) setMethod(generateSpec, B, function(object) {aSpec - object@A@specA; object@specB - 1/aSpec; return(object)}) ## all fine: new(B, value = 2) ## Define class C inheriting from class B ## setClass(C, representation(c = numeric, specC = numeric), contains=c(B)) setMethod(initialize, C, function(.Object, value) {.Object@c - value; .Object - callNextMethod(.Object, value); .Object - generateSpec(.Object); return(.Object)}) setMethod(generateSpec, C, function(object) {bSpec - object@specB; object@specC - 1/bSpec; return(object)}) ## numeric(0) for specB and specC new(C, value = 2) Note, that the generateSpec method of B depends on the specA slot of A, the generateSpec method of class C depends on the specB slot of class B. Note also, that in the initialize method of B one has to give the value argument a default value. When calling the last line of code, specB and specC are numeric(0). As the slot b of the new C class is equal to 2, the initialize method of class B is called but the .Object in the initialize method of class B is an object of class C. And of course in this case method dispatching chooses the generateSpec method of class C - generateSpec for class B gets never called. In general this is exactly what we want, when inheriting classes: dispatching makes OOP programming not only simpler, but actually possible. My question: How can I make this setting work? How can I generate an initialize method for class B that calls ALWAYS the generateSpec method of class B (using as lets me end up in the initialize method of class C with an object of class B as it is only possible to convert a derived object to its base class but not the other way around)? Best Simon __ 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] R number format with Hmisc and knitr
Dear R-Users and R-Devels, I have a problem when using knitr in combination with Hmisc. I generate a data.frame which has mixed scientific and non-scientific numbers inside. In my Latex Table I just want to have non-scientific format, so I call latex(myDataFrame, file = '', cdec = c(0, rep(4, NROW(myDataFrame) - 1)), ) Usually this works, but in this case it doesn't. I do not know why but suspect the mixed data format to be the culprit. What could I do? Using format(, scientific = FALSE) before or options(scipen = 4) before has no influence. Best Simon __ 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] R number format with Hmisc and knitr
Errata: it must say: latex(myDataFrame, file = '', cdec = c(0, rep(4, NCOL(myDataFrame) - 1)) ) But this does not work. Scientific notation is very robust :) Apologize Simon On Jul 31, 2013, at 5:05 PM, Simon Zehnder szehn...@uni-bonn.de wrote: Dear R-Users and R-Devels, I have a problem when using knitr in combination with Hmisc. I generate a data.frame which has mixed scientific and non-scientific numbers inside. In my Latex Table I just want to have non-scientific format, so I call latex(myDataFrame, file = '', cdec = c(0, rep(4, NROW(myDataFrame) - 1)), ) Usually this works, but in this case it doesn't. I do not know why but suspect the mixed data format to be the culprit. What could I do? Using format(, scientific = FALSE) before or options(scipen = 4) before has no influence. Best Simon __ 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-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] MCMClogit: Cannot calculate marginal likelihood with improper prior
Hi, what I see so far is that you have specified your user.prior.density correctly. The error comes from the prior precision matrix B0 in combination with the marginal.likelihood set to Laplace. B0, if not explicitly specified, defaults to zero, which results in eigenvalues of zero. If Laplace is indicated for the marginal.likelihood, the algorithm usually calls an optimization over logpost.logit in BayesianFactors.R where the matrix B0 is tried to be solved by solve(B0) ... as it is a zero matrix its linear equation system is exactly singular and cannot be solved. The Function MCMClogit knows about this fact and gives out a warning Cannot calculate marginal likelihood with improper prior while changing marginal.likelihood to none. So concluding: Choose your user.prior.density with marginal.likelihood = none and all is fine (implicitly it is done so nevertheless). Best Simon P.S. Using a name on a community help list will certainly improve the number of answers to your questions. On Jul 29, 2013, at 3:00 AM, ba0728 haleyb...@att.net wrote: I'm an undergrad who is new to MCMCpack and I haven't been able to find an answer to my problem online yet: I'm attempting to run MCMClogit with a Cauchy proper prior but I'm getting the warning Cannot calculate marginal likelihood with improper prior (my purposes require the marginal likelihood calculation so I understand that I need to use a proper prior). I'm trying to simulate the user-defined independent Cauchy prior with additional args as specified in the MCMCpack User Manual (p. 76, April 2013 version). My input data has been standardized (mean = 0, sd = 0.5 for non-binary variables, and binary variables with mean of 0 and difference of 1 between upper and lower ends) according to the Gelman 2008 paper on logistic regression (www.stat.columbia.edu/~gelman/research/published/priors11.pdf). When I run the example data set (birthwt) from the User Manual, the logpriorfun works correctly allowing the marginal likelihood to be generated. However, when I try running my data with the logprior fun, I get a warning that the prior is improper. Here is the code I am running: *logpriorfun = function(beta, location,scale){ sum(dcauchy(beta, location, scale, log = TRUE)) }* * MCMC.2= MCMClogit(DEAD ~ YEARS + MALE + x1 + x2 + x3+ x4 +x5 + x6 + x7 + x8 + x9, tune= 0.65,burnin =500, mcmc=5000, data = dat, marginal.likelihood = Laplace, user.prior.density=logpriorfun, logfun=TRUE, location = 0, scale=2.5) * *@ The Metropolis acceptance rate was 0.27418 @ Warning message: In MCMClogit(DEAD ~ YEARS + MALE + x1 + x2 + x3 + : Cannot calculate marginal likelihood with improper prior* Any advice on how to fix my arguments so it is a proper prior and will allow me to generate a marginal likelihood using the Laplace approximation? Or how should I be coding a Cauchy proper prior? I'm having problems defining the priors. Thanks, B. -- View this message in context: http://r.789695.n4.nabble.com/MCMClogit-Cannot-calculate-marginal-likelihood-with-improper-prior-tp4672561.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-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] Hmisc ctable rotate option obsolete?
So, I downloaded the source files of Hmisc and changed in the file latex.s line 688 'rotate' to 'sideways'. This does the work for landscape ctables in Latex. I also wrote an email to the package maintainer. I consider this thread as solved. Best Simon On Jul 26, 2013, at 5:34 PM, Simon Zehnder szehn...@uni-bonn.de wrote: Dear R-Users and R-Devels, I may have found a deprecated option for the 'latex' function in the Hmisc package. I am working with Hmisc and knitr and tried the following code: \documentclass{article} \usepackage[utf8]{inputenc} \usepackage{amsmath, amssymb} \usepackage{ctable} %\usepackage{booktabs} \begin{document} results = 'asis'= library(Hmisc) iris.t - head(iris) iris.t[seq(2, NROW(iris.t), by = 2),] - format(iris.t[seq(2, NROW(iris.t), by = 2),], scientific = TRUE) texMat - matrix(, ncol = 5, nrow = 6) texMat[seq(2,nrow(texMat), by = 2), ] - scriptsize latex(iris.t, file = '', landscape = TRUE, dcolumn = TRUE, col.just = c('r','c', 'r','c', 'l'), cdec = c(0, 0, 1, 1, 0), na.blank = TRUE, rowname = '', rowlabel = '', cellTexCmd = texMat, ctable = TRUE, cgroup = c('Observations', ''), n.cgroup = c(4, 1), rgroup = c('',''), n.rgroup = c(3, 3), caption = 'iris' ) @ \end{document} Everything runs fine but the 'landscape' option. It says in the help for 'latex' that if option 'ctable' is set to TRUE the 'rotate' option for ctable is used if 'nadscape' is set TRUE. Looking at the ctable documentary (http://texdoc.net/texmf-dist/doc/latex/ctable/ctable.pdf) in section Change History, I get for version v1.07: General: Added option sideways, option rotate now obsolete. Hasn't this been updated in the Hmisc package? Best Simon __ 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-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] Hmisc ctable rotate option obsolete?
Dear R-Users and R-Devels, I may have found a deprecated option for the 'latex' function in the Hmisc package. I am working with Hmisc and knitr and tried the following code: \documentclass{article} \usepackage[utf8]{inputenc} \usepackage{amsmath, amssymb} \usepackage{ctable} %\usepackage{booktabs} \begin{document} results = 'asis'= library(Hmisc) iris.t - head(iris) iris.t[seq(2, NROW(iris.t), by = 2),] - format(iris.t[seq(2, NROW(iris.t), by = 2),], scientific = TRUE) texMat - matrix(, ncol = 5, nrow = 6) texMat[seq(2,nrow(texMat), by = 2), ] - scriptsize latex(iris.t, file = '', landscape = TRUE, dcolumn = TRUE, col.just = c('r','c', 'r','c', 'l'), cdec = c(0, 0, 1, 1, 0), na.blank = TRUE, rowname = '', rowlabel = '', cellTexCmd = texMat, ctable = TRUE, cgroup = c('Observations', ''), n.cgroup = c(4, 1), rgroup = c('',''), n.rgroup = c(3, 3), caption = 'iris' ) @ \end{document} Everything runs fine but the 'landscape' option. It says in the help for 'latex' that if option 'ctable' is set to TRUE the 'rotate' option for ctable is used if 'nadscape' is set TRUE. Looking at the ctable documentary (http://texdoc.net/texmf-dist/doc/latex/ctable/ctable.pdf) in section Change History, I get for version v1.07: General: Added option sideways, option rotate now obsolete. Hasn't this been updated in the Hmisc package? Best Simon __ 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] Check the class of an object
Dear R-Users and R-Devels, I have large project based on S4 classes. While writing my unit tests I found out, that 'is' cannot test for a specific class, as also inherited classes can be treated as their super classes. I need to do checks for specific classes. What I do right now is sth. like if (class(myClass) == firstClass) { } else if (class(myClass) == secondClass) { } Is this the usual way how classes are checked in R? I was expecting some specific method (and 'inherits' or 'extends' is not what I look for)... Best Simon [[alternative HTML version deleted]] __ 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] Check the class of an object
Hi David, thanks for the reply. You are right. Using the %in% is more stable and I gonna change my code. When testing for a specific class using 'is' one has to start at the lowest heir and walk up the inheritance structure. Starting at the checks at the root will always give TRUE. Having a structure which is quite complicated let me move to the check I suggested in my first mail. Best Simon On Jul 23, 2013, at 6:15 PM, David Winsemius dwinsem...@comcast.net wrote: On Jul 23, 2013, at 5:36 AM, Simon Zehnder wrote: Dear R-Users and R-Devels, I have large project based on S4 classes. While writing my unit tests I found out, that 'is' cannot test for a specific class, as also inherited classes can be treated as their super classes. I need to do checks for specific classes. What I do right now is sth. like if (class(myClass) == firstClass) { I would think that you would need to use `%in%` instead. if( firstClass %in% class(myObject) ){ Objects can have more than one class, so testing with == would fail in those instances. } else if (class(myClass) == secondClass) { } Is this the usual way how classes are checked in R? Well, `inherits` IS the usual way. I was expecting some specific method (and 'inherits' or 'extends' is not what I look for)... Best Simon [[alternative HTML version deleted]] Plain-text format is the recommended format for Rhelp -- David Winsemius Alameda, CA, USA __ 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] Check the class of an object
Hi Martin, I didn't know that. But that is even more comfortable for checking. Thanks for the quick update! Best Simon P.S. And thanks for the online documents about S4 - I could already learn a lot! On Jul 23, 2013, at 7:11 PM, Martin Morgan mtmor...@fhcrc.org wrote: On 07/23/2013 09:59 AM, Simon Zehnder wrote: Hi David, thanks for the reply. You are right. Using the %in% is more stable and I gonna change my code. you said you were you were using S4 classes. S4 classes do not report vectors of length != 1, from ?class For objects which have a formal class, its name is returned by 'class' as a character vector of length one so a first unit test could be stopifnot(length(class(myObject)) != 1L) When testing for a specific class using 'is' one has to start at the lowest heir and walk up the inheritance structure. Starting at the checks at the root will always give TRUE. Having a structure which is quite complicated let me move to the check I suggested in my first mail. Best Simon On Jul 23, 2013, at 6:15 PM, David Winsemius dwinsem...@comcast.net wrote: On Jul 23, 2013, at 5:36 AM, Simon Zehnder wrote: Dear R-Users and R-Devels, I have large project based on S4 classes. While writing my unit tests I found out, that 'is' cannot test for a specific class, as also inherited classes can be treated as their super classes. I need to do checks for specific classes. What I do right now is sth. like if (class(myClass) == firstClass) { I would think that you would need to use `%in%` instead. if( firstClass %in% class(myObject) ){ Objects can have more than one class, so testing with == would fail in those instances. } else if (class(myClass) == secondClass) { } Is this the usual way how classes are checked in R? Well, `inherits` IS the usual way. I was expecting some specific method (and 'inherits' or 'extends' is not what I look for)... Best Simon [[alternative HTML version deleted]] Plain-text format is the recommended format for Rhelp -- David Winsemius Alameda, CA, USA __ 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. -- Computational Biology / Fred Hutchinson Cancer Research Center 1100 Fairview Ave. N. PO Box 19024 Seattle, WA 98109 Location: Arnold Building M1 B861 Phone: (206) 667-2793 __ 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] Check the class of an object
Hi Hervé, thank you very much for your reply! This makes the different treatment of S3 and S4 objects by 'class' clear. Best Simon On Jul 23, 2013, at 7:20 PM, Hervé Pagès hpa...@fhcrc.org wrote: Hi, On 07/23/2013 09:59 AM, Simon Zehnder wrote: Hi David, thanks for the reply. You are right. Using the %in% is more stable and I gonna change my code. Unlike with S3 objects, class() on an S4 object can only return 1 class. Also note that, on an S3 object, doing firstClass %in% class(myObject) is equivalent to doing inherits(myObject, firstClass), which is what you said you wanted to avoid. The most specific class should be the first so if that's what you wanted to check, you could do class(myObject)[1] == firstClass But that precaution is not needed if 'myObject' is guaranteed to be an S4 object (although when writing a unit test, one should probably discard any guarantee of that sort). Cheers, H. When testing for a specific class using 'is' one has to start at the lowest heir and walk up the inheritance structure. Starting at the checks at the root will always give TRUE. Having a structure which is quite complicated let me move to the check I suggested in my first mail. Best Simon On Jul 23, 2013, at 6:15 PM, David Winsemius dwinsem...@comcast.net wrote: On Jul 23, 2013, at 5:36 AM, Simon Zehnder wrote: Dear R-Users and R-Devels, I have large project based on S4 classes. While writing my unit tests I found out, that 'is' cannot test for a specific class, as also inherited classes can be treated as their super classes. I need to do checks for specific classes. What I do right now is sth. like if (class(myClass) == firstClass) { I would think that you would need to use `%in%` instead. if( firstClass %in% class(myObject) ){ Objects can have more than one class, so testing with == would fail in those instances. } else if (class(myClass) == secondClass) { } Is this the usual way how classes are checked in R? Well, `inherits` IS the usual way. I was expecting some specific method (and 'inherits' or 'extends' is not what I look for)... Best Simon [[alternative HTML version deleted]] Plain-text format is the recommended format for Rhelp -- David Winsemius Alameda, CA, USA __ 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. -- Hervé Pagès Program in Computational Biology Division of Public Health Sciences Fred Hutchinson Cancer Research Center 1100 Fairview Ave. N, M1-B514 P.O. Box 19024 Seattle, WA 98109-1024 E-mail: hpa...@fhcrc.org Phone: (206) 667-5791 Fax:(206) 667-1319 __ 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] S4 method signature - integer matrix
Hi David, thanks for the reply! The prototype solution is I think not the appropriate solution, when the method adds the integer matrix next to the class object into the argument list. Also the .validObject method can only be applied to the object for which the method is called or I have to write an integer matrix class that has its own .validObject method. From your suggestions I had a direction to search more precisely and I found this: Set the signature to matrix: checks for matrix when method is called. Then check for typeof(mymatrix) == integer: Return an exception if FALSE. In case of an index I can add then an all(apply(mymatrix, c(1, 2), , 0)): Return an exception if FALSE Thank you very much for your advice. It helped me to find a good solution. Best Simon On Jul 20, 2013, at 12:03 AM, David Winsemius dwinsem...@comcast.net wrote: On Jul 19, 2013, at 9:54 AM, Simon Zehnder wrote: Dear R-Users and R-Devels, I am programming a package with S4 classes and I search for a solution of the following problem: If you want an S4 method to await an integer argument you set the signature like this setMethod(myfunction, signature(object = myClass, y = integer), function(object, y) {//Do sth}) Now, if you want the method to await an integer matrix or array there is only one way how to define it setMethod(myfunction, signature(object = myClass, y = matrix), function(object, y) {//Do sth}) or setMethod(myfunction, signature(object = myClass, y = array), function(object, y) {//Do sth}) am I right? I don't think that is the only way. You could also test for mode being 'numeric' at the signature level and then within the validation check to see whether it has dimensions. See ?is.numeric and ?validObject There is also the option of using a prototype. As I understand the S4 checks they will apply an is.mode or is.class test as long as there is an correspond function that returns a logical. see ... ?is WIth this you can also pass a numeric matrix. How would you check for an integer matrix in an S4 method (in the index function of R I think it is just coerced to integer)? It would need be both integer mode... ?is.integer and have dimensions... ?dim Furthermore: How would you check for an integer matrix with values bigger one (so the typical R indices)? Just add a test in the 'validity' code. See ... ?setClass Is there a way how it is usually done in R (I think probably with apply())? Or is it usually better to throw an exception? I don't see a reason to use apply. -- David Winsemius Alameda, CA, USA __ 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] S4 method signature - integer matrix
Dear R-Users and R-Devels, I am programming a package with S4 classes and I search for a solution of the following problem: If you want an S4 method to await an integer argument you set the signature like this setMethod(myfunction, signature(object = myClass, y = integer), function(object, y) {//Do sth}) Now, if you want the method to await an integer matrix or array there is only one way how to define it setMethod(myfunction, signature(object = myClass, y = matrix), function(object, y) {//Do sth}) or setMethod(myfunction, signature(object = myClass, y = array), function(object, y) {//Do sth}) am I right? WIth this you can also pass a numeric matrix. How would you check for an integer matrix in an S4 method (in the index function of R I think it is just coerced to integer)? Furthermore: How would you check for an integer matrix with values bigger one (so the typical R indices)? Is there a way how it is usually done in R (I think probably with apply())? Or is it usually better to throw an exception? Best Simon __ 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] 'save' method for S4 class
Hi Christopher, I think, that save is no generic function like plot, show, etc. So at first you have to determine a generic. setGeneric(save, function(x, file_Path) standardGeneric(save)) Now your definition via setMethod. Best Simon On Jul 18, 2013, at 12:09 PM, Christofer Bogaso bogaso.christo...@gmail.com wrote: Hello again, I am trying to define the 'save' method for my S4 class as below: setClass(MyClass, representation( Slot1 = data.frame )) setMethod(save, MyClass, definition = function(x, file_Path) { write.table(x@Slot1, file = file_Path, append = FALSE, quote = TRUE, sep = ,, eol = \n, na = NA, dec = ., row.names = FALSE, col.names = TRUE, qmethod = c(escape, double), fileEncoding = ) }) However while doing this I am getting following error: Error in conformMethod(signature, mnames, fnames, f, fdef, definition) : in method for ‘save’ with signature ‘list=MyClass’: formal arguments (list = MyClass, file = MyClass, ascii = MyClass, version = MyClass, envir = MyClass, compress = MyClass, compression_level = MyClass, eval.promises = MyClass, precheck = MyClass) omitted in the method definition cannot be in the signature Can somebody point me what will be the correct approach to define 'save' method for S4 class? Thanks and regards, __ 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-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] 'save' method for S4 class
Hi Christofer, R's S4 class system builds on generics, that is one class defines the raw frame of a function and maybe an own distinct implementation and other classes, that (usually) inherit from this class, set their special implementations of it. Consider a generic function in R as being an abstract function, just telling inheriting classes how to build such a function, i.e. what to put in. Furthermore it enables the dispatches like callNextMethod, which calls the method from the parent (be aware if you have more than one direct parent!). Under the hood, R must know from the input arguments what function should be called. If you put in an B class that inherits from an A class it must know, what to do with a user-defined save function when the input argument is for example an object of class B. If you put in an A class object, it should know how to proceed this object and so on. So the generic method gives a kind of internal map for R how to call this function in which contexts. The standardGeneric() function: We must provide a name and a definition of the generic function. In most cases you just want to tell, what arguments the function should expect. With the function standardGeneric(name) R creates a generic which can then be dispatched in form of distinct methods. So if you call save R sees, that it is a generic function and checks which of the distinct functions to call, evaluating the arguments. If it finds no specific function for the arguments at hand, it calls a default function. If, as in your case, the function name is already taken, R sets, when you use setGeneric, the usual save function as the default method, getting called when nothing else fits ... resulting in either an error because of mismatching function arguments or in saving in a way you do not want to. Solution is here to define inside the setGeneric also the default function via the argument useAsDefault = function() {}. So generic function (setGeneric): What should be done? Distinct function (setMethod): How should it be done? You find a very good introduction to the S4 class system here: http://cran.r-project.org/doc/contrib/Genolini-S4tutorialV0-5en.pdf Deep down information can be found for example in the books: Programming with Data and Software for Data Analysis and Computing with R, both from John Chambers. Best Simon On Jul 18, 2013, at 8:23 PM, Christofer Bogaso bogaso.christo...@gmail.com wrote: Hi Simon, Thanks for your pointer. However could you please explain what 'function(x, file_Path) standardGeneric(save)' means here? The underlying help files look quite rocket science for me! Thanks for your time. Thanks and regards, On Thu, Jul 18, 2013 at 4:32 PM, Simon Zehnder szehn...@uni-bonn.de wrote: Hi Christopher, I think, that save is no generic function like plot, show, etc. So at first you have to determine a generic. setGeneric(save, function(x, file_Path) standardGeneric(save)) Now your definition via setMethod. Best Simon On Jul 18, 2013, at 12:09 PM, Christofer Bogaso bogaso.christo...@gmail.com wrote: Hello again, I am trying to define the 'save' method for my S4 class as below: setClass(MyClass, representation( Slot1 = data.frame )) setMethod(save, MyClass, definition = function(x, file_Path) { write.table(x@Slot1, file = file_Path, append = FALSE, quote = TRUE, sep = ,, eol = \n, na = NA, dec = ., row.names = FALSE, col.names = TRUE, qmethod = c(escape, double), fileEncoding = ) }) However while doing this I am getting following error: Error in conformMethod(signature, mnames, fnames, f, fdef, definition) : in method for save with signature list=MyClass: formal arguments (list = MyClass, file = MyClass, ascii = MyClass, version = MyClass, envir = MyClass, compress = MyClass, compression_level = MyClass, eval.promises = MyClass, precheck = MyClass) omitted in the method definition cannot be in the signature Can somebody point me what will be the correct approach to define 'save' method for S4 class? Thanks and regards, __ 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. [[alternative HTML version deleted]] __ 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] Looking for knitr example for beginner (NO RStudio)
Hi Mike, I found my way with this little blog: http://yihui.name/knitr/demo/editors/ The .Rnw files are created very well in a Latex editor. Everything else can be easily googled. The command via knitr::knit2pdf works very fine if you use the chunks. If you are trying to compile an Rtex file, this I do not know either (I like the symbols though in for example https://github.com/yihui/knitr-examples/blob/master/005-latex.Rtex). But the .Rnw files are compiled pretty nice in e.g. texmaker, as described in the blog. Use for example this source file: https://github.com/yihui/knitr/blob/master/inst/examples/knitr-minimal.Rnw Hope this helps Best Simon On Jul 18, 2013, at 8:52 PM, C W tmrs...@gmail.com wrote: How do you create a .Rnw file, in R or LaTex? I don't think any tutorial mentions it. btw, I am very new to the terms like markdown, so I don't understand markdown to HTML. I am reading here http://biostat.mc.vanderbilt.edu/wiki/Main/KnitrHowto that you need to compile at terminal. I do not know terminal, is there other ways? Could you do a video on just simple R? I have seen 3 videos on R Studio including yours. Mike On Thu, Jul 18, 2013 at 2:43 PM, Yihui Xie x...@yihui.name wrote: I'm not sure what your question really is. You do not have to use RStudio, but it will be much easier to get started with RStudio, because it does a lot of automatic conversion behind the scenes (e.g. tex to PDF, markdown to HTML, ...). If you want a pure solution without any text editor support, the answer is library(knitr) knit('your_input_file') For example, knit('foo.Rnw') gives you foo.tex; if you are familiar with LaTeX, you can mess with this foo.tex now (outside of R). Minimal examples for different document formats are at http://yihui.name/knitr/demo/minimal/ (you must have read this page), and more examples at https://github.com/yihui/knitr-examples If you are asking about the internals of knitr, Luke, use the source: https://github.com/yihui/knitr Or for a more comprehensive introduction, see http://www.crcpress.com/product/isbn/9781482203530 Regards, Yihui -- Yihui Xie xieyi...@gmail.com Phone: 206-667-4385 Web: http://yihui.name Fred Hutchinson Cancer Research Center, Seattle On Thu, Jul 18, 2013 at 11:13 AM, C W tmrs...@gmail.com wrote: Hi everyone, I am using package knitr, FIRST TIME. I don't have access to RStudio. Read through Yihui's page, didn't find it helpful. Stuck on terms Rnw, GFM (GitHub Flavored Markdown). Never used Sweave, so the reference is not helping. Is there a simple step-by-step example WITHOUT RStudio? My question: What is the procedure? The documentation explains the functions, but does not say how to operate between R and LaTex. Mike __ 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-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-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] Looking for knitr example for beginner (NO RStudio)
The executable is in case of knitr always Rscript. On a mac it is simply Rscript on windows it is Rscript.exe. This should be on your PATH. If you are not sure, open the Mac Terminal and type Rscript --version. If it does not say Command not found all is fine. Best Simon On Jul 18, 2013, at 10:09 PM, C W tmrs...@gmail.com wrote: http://tex.stackexchange.com/questions/85154/knitr-with-texworks/85165#85165 In step 3: add the executable file (step 3). What is the executable file? Locate package knitr directory path in R? Mike On Thu, Jul 18, 2013 at 3:56 PM, C W tmrs...@gmail.com wrote: Actually, I see it at the bottom. Sorry! On Thu, Jul 18, 2013 at 3:44 PM, C W tmrs...@gmail.com wrote: Hi Simon, I am on OS X Lion, I have TeXworks, I don't have knitr as an option. How do I install that into TeXworks? Seems like I have to something in terminal? Mike On Thu, Jul 18, 2013 at 3:31 PM, Simon Zehnder szehn...@uni-bonn.de wrote: Hi Mike, I found my way with this little blog: http://yihui.name/knitr/demo/editors/ The .Rnw files are created very well in a Latex editor. Everything else can be easily googled. The command via knitr::knit2pdf works very fine if you use the chunks. If you are trying to compile an Rtex file, this I do not know either (I like the symbols though in for example https://github.com/yihui/knitr-examples/blob/master/005-latex.Rtex). But the .Rnw files are compiled pretty nice in e.g. texmaker, as described in the blog. Use for example this source file: https://github.com/yihui/knitr/blob/master/inst/examples/knitr-minimal.Rnw Hope this helps Best Simon On Jul 18, 2013, at 8:52 PM, C W tmrs...@gmail.com wrote: How do you create a .Rnw file, in R or LaTex? I don't think any tutorial mentions it. btw, I am very new to the terms like markdown, so I don't understand markdown to HTML. I am reading here http://biostat.mc.vanderbilt.edu/wiki/Main/KnitrHowto that you need to compile at terminal. I do not know terminal, is there other ways? Could you do a video on just simple R? I have seen 3 videos on R Studio including yours. Mike On Thu, Jul 18, 2013 at 2:43 PM, Yihui Xie x...@yihui.name wrote: I'm not sure what your question really is. You do not have to use RStudio, but it will be much easier to get started with RStudio, because it does a lot of automatic conversion behind the scenes (e.g. tex to PDF, markdown to HTML, ...). If you want a pure solution without any text editor support, the answer is library(knitr) knit('your_input_file') For example, knit('foo.Rnw') gives you foo.tex; if you are familiar with LaTeX, you can mess with this foo.tex now (outside of R). Minimal examples for different document formats are at http://yihui.name/knitr/demo/minimal/ (you must have read this page), and more examples at https://github.com/yihui/knitr-examples If you are asking about the internals of knitr, Luke, use the source: https://github.com/yihui/knitr Or for a more comprehensive introduction, see http://www.crcpress.com/product/isbn/9781482203530 Regards, Yihui -- Yihui Xie xieyi...@gmail.com Phone: 206-667-4385 Web: http://yihui.name Fred Hutchinson Cancer Research Center, Seattle On Thu, Jul 18, 2013 at 11:13 AM, C W tmrs...@gmail.com wrote: Hi everyone, I am using package knitr, FIRST TIME. I don't have access to RStudio. Read through Yihui's page, didn't find it helpful. Stuck on terms Rnw, GFM (GitHub Flavored Markdown). Never used Sweave, so the reference is not helping. Is there a simple step-by-step example WITHOUT RStudio? My question: What is the procedure? The documentation explains the functions, but does not say how to operate between R and LaTex. Mike __ 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-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-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] Looking for knitr example for beginner (NO RStudio)
Hi Mike, if you browse the folders, you find always the Rscript binary (the executable) under /Library/Frameworks/R.framework/Versions/.../Resources/Rscript. Do not forget to give your tex file the extension .Rnw! Then surround each Rcode with write here a nameadd later further options (important one; results = 'asis')= Here your r code as you do it in the R shell at the end a @. Always inside the \begin{document} \end{document} tags. Best Simon On Jul 18, 2013, at 10:22 PM, C W tmrs...@gmail.com wrote: Thanks, Simon. I would never figured it out! I apologize if I sound frustrated, because I am. @package author: you have a great package, but I think a lot of the directions are hand waving. For the newbies, this leads to more confusion. @Berend: I am using OS X. Mike On Thu, Jul 18, 2013 at 4:19 PM, Berend Hasselman b...@xs4all.nl wrote: On 18-07-2013, at 22:09, C W tmrs...@gmail.com wrote: http://tex.stackexchange.com/questions/85154/knitr-with-texworks/85165#85165 In step 3: add the executable file (step 3). What is the executable file? Locate package knitr directory path in R? From the window: Executable == Program. So the executable is Rscript.exe. Berend __ 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] Serialize data.frame to database
Hi Jeff, I think you are more right than me. I have not much experience with R-specific objects and Databases. I just know, that for languages like Java, C, etc. this process is a usual one for storing for example matrices (in case of course that I do not have to access specific elements inside of the matrix but only the matrix as a whole). Thank you for the tip with the R-sig-DB list. I switch over to this list. Best Simon On Jul 16, 2013, at 2:34 AM, Jeff Newmiller jdnew...@dcn.davis.ca.us wrote: I could be wrong, but I would guess that doing what you are describing is very unusual. Most of the time the data frame is mapped to a table in the database so the rows can be searched. Storing data frames as BLOBs really seems odd. Note that there is an R-sig-db mailing list for questions of this type. --- Jeff NewmillerThe . . Go Live... DCN:jdnew...@dcn.davis.ca.usBasics: ##.#. ##.#. Live Go... Live: OO#.. Dead: OO#.. Playing Research Engineer (Solar/BatteriesO.O#. #.O#. with /Software/Embedded Controllers) .OO#. .OO#. rocks...1k --- Sent from my phone. Please excuse my brevity. Simon Zehnder szehn...@uni-bonn.de wrote: Dear R-Users, I need a very fast and reliable database solution so I try to serialize a data.frame (to binary data) and to store this data to an SQLite database. This is what I tried to do: library(RSQLite) con - dbDriver(SQLite) db - dbConnect(con, test) dbSendQuery(db, 'CREATE TABLE frames(simID INT, data BLOB)') data.bin - serialize(iris, NULL, ascii = FALSE) dbSendQuery(db, paste(INSERT INTO frames VALUES(1, X', data.bin, '), sep = )) data.bin2 - dbGetQuery(db, SELECT DATA FROM frames WHERE simID = 1) data.bin2 data 1 58 So, only the first entry of data.bin is saved to the database. I tried to first convert the binary data to raw data: data.raw - rawToChar(data.bin) Error in rawToChar(data.bin) : embedded nul in string: 'X\n\0\0\0\002\0\003\0\001\0\002\003\0\0\0\003\023\0\0\0\005\0\0\0\016\0\0\0\x96@\024ff@\023\x99\x99\x99\x99\x99\x9a@\022\xcc\xcc\xcc\xcc\xcc\xcd@\022ff@\024\0\0\0\0\0\0@\025\x99\x99\x99\x99\x99\x9a@\022ff@\024\0\0\0\0\0\0@\021\x99\x99\x99\x99\x99\x9a@\023\x99\x99\x99\x99\x99\x9a@\025\x99\x99\x99\x99\x99\x9a@\02333@\02333@\02133@\02733@\026\xcc\xcc\xcc\xcc\xcc\xcd@\025\x99\x99\x99\x99\x99\x9a@\024ff@\026\xcc\xcc\xcc\xcc\xcc\xcd@\024ff@\025\x99\x99\x99\x99\x99\x9a@\024ff@\022ff@\024ff@\02333@\024\0\0\0\0\0\0@\024\0\0\0\0\0\0@\024\xcc\xcc\xcc\xcc\xcc\xcd@\024\xcc\xcc\xcc\xcc\xcc\xcd@\022\xcc\xcc\xcc\xcc\xcc\xcd@\02333@\025\x99\x99\x99\x99\x99\x9a@\024\xcc\xcc\xcc\xcc\xcc\xcd@\026\0\0\0\0\0\0@\023\x99\x99\x99\x99\x99\x9a@\024\0\0\0\0\0\0@\026\0\0\0\0\0\0@\023\x99\x99\x99\x99\x99\x9a@\021\x99\x99\x99\x99\x99\x9a@\024ff@\024\0\0\0\0\0\0@\022\0\0\0\0\0\0@\021\x99\x99\x99\x99\x99\x9a@\024\0\0\0\0\! 0\0 I don't know what this error should tell me. Then I tried to use the ASCII format data.ascii - serialize(iris, NULL, ascii = TRUE) data.raw - rawToChar(data.ascii) dbSendQuery(db, DELETE FROM frames) dbSendQuery(db, paste(INSERT INTO frames VALUES(1, X', data.raw, '), sep = )) Error in sqliteExecStatement(conn, statement, ...) : RS-DBI driver: (error in statement: unrecognized token: X'A This also does not work. It seems the driver does not deal that nicely with the regular INSERT query for BLOB objects in SQLite. Then I used a simpler way: dbSendQuery(db, DELETE FROM frames) dbSendQuery(db, DROP TABLE frames) dbSendQuery(db, 'CREATE TABLE frames(simID INT, data TEXT DEFAULT NULL)') dbSendQuery(db, paste(INSERT INTO frames VALUES(1, ', data.raw, '), sep = )) data.bin2 - dbGetQuery(db, SELECT data FROM frames WHERE simID = 1) Nice, that worked. Now I want to unserialize the data: unserialize(data.bin2) Error in unserialize(data.bin2) : 'connection' must be a connection unserialize(data.bin2[1, 'data']) Error in unserialize(data.bin2[1, data]) : character vectors are no longer accepted by unserialize() I feel a little stuck here, but I am very sure, that converting data.frames to binary data and storing them to a database is not that unusual. So I hope somebody has already done this and could give me the missing piece. Best Simon __ 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-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r
Re: [R] Serialize data.frame to database
Hi Rainer, dbWriteTable is a nice function but in my case I need something that can actually save a dataframe in one row of a table. That is why I want to serialize my data.frame. Best Simon On Jul 16, 2013, at 3:05 PM, Rainer Schuermann rainer.schuerm...@gmx.net wrote: Maybe a simple dbWriteTable( db, frames, iris ) does what you want? On Monday 15 July 2013 23:43:18 Simon Zehnder wrote: Dear R-Users, I need a very fast and reliable database solution so I try to serialize a data.frame (to binary data) and to store this data to an SQLite database. This is what I tried to do: library(RSQLite) con - dbDriver(SQLite) db - dbConnect(con, test) dbSendQuery(db, 'CREATE TABLE frames(simID INT, data BLOB)') data.bin - serialize(iris, NULL, ascii = FALSE) dbSendQuery(db, paste(INSERT INTO frames VALUES(1, X', data.bin, '), sep = )) data.bin2 - dbGetQuery(db, SELECT DATA FROM frames WHERE simID = 1) data.bin2 data 1 58 So, only the first entry of data.bin is saved to the database. I tried to first convert the binary data to raw data: data.raw - rawToChar(data.bin) Error in rawToChar(data.bin) : embedded nul in string: 'X\n\0\0\0\002\0\003\0\001\0\002\003\0\0\0\003\023\0\0\0\005\0\0\0\016\0\0\0\x96@\024ff@\023\x99\x99\x99\x99\x99\x9a@\022\xcc\xcc\xcc\xcc\xcc\xcd@\022ff@\024\0\0\0\0\0\0@\025\x99\x99\x99\x99\x99\x9a@\022ff@\024\0\0\0\0\0\0@\021\x99\x99\x99\x99\x99\x9a@\023\x99\x99\x99\x99\x99\x9a@\025\x99\x99\x99\x99\x99\x9a@\02333@\02333@\02133@\02733@\026\xcc\xcc\xcc\xcc\xcc\xcd@\025\x99\x99\x99\x99\x99\x9a@\024ff@\026\xcc\xcc\xcc\xcc\xcc\xcd@\024ff@\025\x99\x99\x99\x99\x99\x9a@\024ff@\022ff@\024ff@\02333@\024\0\0\0\0\0\0@\024\0\0\0\0\0\0@\024\xcc\xcc\xcc\xcc\xcc\xcd@\024\xcc\xcc\xcc\xcc\xcc\xcd@\022\xcc\xcc\xcc\xcc\xcc\xcd@\02333@\025\x99\x99\x99\x99\x99\x9a@\024\xcc\xcc\xcc\xcc\xcc\xcd@\026\0\0\0\0\0\0@\023\x99\x99\x99\x99\x99\x9a@\024\0\0\0\0\0\0@\026\0\0\0\0\0\0@\023\x99\x99\x99\x99\x99\x9a@\021\x99\x99\x99\x99\x99\x9a@\024ff@\024\0\0\0\0\0\0@\022\0\0\0\0\0\0@\021\x99\x99\x99\x99\x99\x9a@\024\0\0\0\! 0\! 0\0 I don't know what this error should tell me. Then I tried to use the ASCII format data.ascii - serialize(iris, NULL, ascii = TRUE) data.raw - rawToChar(data.ascii) dbSendQuery(db, DELETE FROM frames) dbSendQuery(db, paste(INSERT INTO frames VALUES(1, X', data.raw, '), sep = )) Error in sqliteExecStatement(conn, statement, ...) : RS-DBI driver: (error in statement: unrecognized token: X'A This also does not work. It seems the driver does not deal that nicely with the regular INSERT query for BLOB objects in SQLite. Then I used a simpler way: dbSendQuery(db, DELETE FROM frames) dbSendQuery(db, DROP TABLE frames) dbSendQuery(db, 'CREATE TABLE frames(simID INT, data TEXT DEFAULT NULL)') dbSendQuery(db, paste(INSERT INTO frames VALUES(1, ', data.raw, '), sep = )) data.bin2 - dbGetQuery(db, SELECT data FROM frames WHERE simID = 1) Nice, that worked. Now I want to unserialize the data: unserialize(data.bin2) Error in unserialize(data.bin2) : 'connection' must be a connection unserialize(data.bin2[1, 'data']) Error in unserialize(data.bin2[1, data]) : character vectors are no longer accepted by unserialize() I feel a little stuck here, but I am very sure, that converting data.frames to binary data and storing them to a database is not that unusual. So I hope somebody has already done this and could give me the missing piece. Best Simon __ 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. - - - - - Der NSA keine Chance: e-mail verschluesseln! http://www.gpg4win.org/ __ 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-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] Setting Derived Class Slots
Hi Steve, it seems to me, that you want to pass an object inside of getDelayProfile. In this case you must use setReplaceMethod(getDelayProfile-,) in your definition inside the virtual and of course also in the specification of OP_Appt. R does not know, that your function should give back an object, that replaces the one on which it was called. Otherwise making any modifications on the slot of the object inside the function will never arrive at the object outside of the function. They have two different environments and as soon as you modify an object inside of a function R builds a copy of it. See more about this by googling R pass-by-value. By the way: I wouldn't call this function getDelayProfile as getters and setters have usually differing a different meaning (encapsulation) in OO programming. Use rather something like fetchDelayProfile or initializeDelayProfile. It will probably make it easier for Users and Programmers to work with your program. Hope this helps Best Simon On Jul 16, 2013, at 3:36 PM, Steve Creamer stephen.crea...@nhs.net wrote: [[alternative HTML version deleted]] __ 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] Fitting Mixture distributions
Hi Tjun Kiat Teo, you try to fit a Normal mixture to some data. The Normal mixture is very delicate when it comes to parameter search: If the variance gets closer and closer to zero, the log Likelihood becomes larger and larger for any values of the remaining parameters. Furthermore for the EM algorithm it is known, that it takes sometimes very long until convergence is reached. Try the following: Use as starting values for the component parameters: start.par - mean(your.data, na.rm = TRUE) + sd(your.data, na.rm = TRUE) * runif(K) For the weights just use either 1/K or the R cluster function with K clusters Here K is the number of components. Further enlarge the maximum number of iterations. What you could also try is to randomize start parameters and run an SEM (Stochastic EM). In my opinion the better method is in this case a Bayesian method: MCMC. Best Simon On Jul 16, 2013, at 10:59 PM, Tjun Kiat Teo teotj...@gmail.com wrote: I was trying to use the normixEM in mixtools and I got this error message. And I got this error message One of the variances is going to zero; trying new starting values. Error in normalmixEM(as.matrix(temp[[gc]][, -(f + 1)])) : Too many tries! Are there any other packages for fitting mixture distributions ? Tjun Kiat Teo [[alternative HTML version deleted]] __ 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-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] Serialize data.frame to database
Dear R-Users, I need a very fast and reliable database solution so I try to serialize a data.frame (to binary data) and to store this data to an SQLite database. This is what I tried to do: library(RSQLite) con - dbDriver(SQLite) db - dbConnect(con, test) dbSendQuery(db, 'CREATE TABLE frames(simID INT, data BLOB)') data.bin - serialize(iris, NULL, ascii = FALSE) dbSendQuery(db, paste(INSERT INTO frames VALUES(1, X', data.bin, '), sep = )) data.bin2 - dbGetQuery(db, SELECT DATA FROM frames WHERE simID = 1) data.bin2 data 1 58 So, only the first entry of data.bin is saved to the database. I tried to first convert the binary data to raw data: data.raw - rawToChar(data.bin) Error in rawToChar(data.bin) : embedded nul in string: 'X\n\0\0\0\002\0\003\0\001\0\002\003\0\0\0\003\023\0\0\0\005\0\0\0\016\0\0\0\x96@\024ff@\023\x99\x99\x99\x99\x99\x9a@\022\xcc\xcc\xcc\xcc\xcc\xcd@\022ff@\024\0\0\0\0\0\0@\025\x99\x99\x99\x99\x99\x9a@\022ff@\024\0\0\0\0\0\0@\021\x99\x99\x99\x99\x99\x9a@\023\x99\x99\x99\x99\x99\x9a@\025\x99\x99\x99\x99\x99\x9a@\02333@\02333@\02133@\02733@\026\xcc\xcc\xcc\xcc\xcc\xcd@\025\x99\x99\x99\x99\x99\x9a@\024ff@\026\xcc\xcc\xcc\xcc\xcc\xcd@\024ff@\025\x99\x99\x99\x99\x99\x9a@\024ff@\022ff@\024ff@\02333@\024\0\0\0\0\0\0@\024\0\0\0\0\0\0@\024\xcc\xcc\xcc\xcc\xcc\xcd@\024\xcc\xcc\xcc\xcc\xcc\xcd@\022\xcc\xcc\xcc\xcc\xcc\xcd@\02333@\025\x99\x99\x99\x99\x99\x9a@\024\xcc\xcc\xcc\xcc\xcc\xcd@\026\0\0\0\0\0\0@\023\x99\x99\x99\x99\x99\x9a@\024\0\0\0\0\0\0@\026\0\0\0\0\0\0@\023\x99\x99\x99\x99\x99\x9a@\021\x99\x99\x99\x99\x99\x9a@\024ff@\024\0\0\0\0\0\0@\022\0\0\0\0\0\0@\021\x99\x99\x99\x99\x99\x9a@\024\0\0\0\0\! 0\0 I don't know what this error should tell me. Then I tried to use the ASCII format data.ascii - serialize(iris, NULL, ascii = TRUE) data.raw - rawToChar(data.ascii) dbSendQuery(db, DELETE FROM frames) dbSendQuery(db, paste(INSERT INTO frames VALUES(1, X', data.raw, '), sep = )) Error in sqliteExecStatement(conn, statement, ...) : RS-DBI driver: (error in statement: unrecognized token: X'A This also does not work. It seems the driver does not deal that nicely with the regular INSERT query for BLOB objects in SQLite. Then I used a simpler way: dbSendQuery(db, DELETE FROM frames) dbSendQuery(db, DROP TABLE frames) dbSendQuery(db, 'CREATE TABLE frames(simID INT, data TEXT DEFAULT NULL)') dbSendQuery(db, paste(INSERT INTO frames VALUES(1, ', data.raw, '), sep = )) data.bin2 - dbGetQuery(db, SELECT data FROM frames WHERE simID = 1) Nice, that worked. Now I want to unserialize the data: unserialize(data.bin2) Error in unserialize(data.bin2) : 'connection' must be a connection unserialize(data.bin2[1, 'data']) Error in unserialize(data.bin2[1, data]) : character vectors are no longer accepted by unserialize() I feel a little stuck here, but I am very sure, that converting data.frames to binary data and storing them to a database is not that unusual. So I hope somebody has already done this and could give me the missing piece. Best Simon __ 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] Set window title for plot on any OS
Haha, good point! So, the correct code is: any(names(formals(getOption(device))) == title) { dev.new(title = title) } Thanks for correcting my approach Rolf! Best Simon On Jul 14, 2013, at 1:02 AM, Rolf Turner rolf.tur...@xtra.co.nz wrote: I think you ought to take a look at fortune(Lewis Carroll) cheers, Rolf Turner On 14/07/13 01:45, Simon Zehnder wrote: Hi Duncan, thank you very much for your advice! That makes it all work. I check in addition for a title argument in the device via if (any(names(getOption(device)) == TRUE)) { dev.new(title = title) } Thanks again! Simon On Jul 13, 2013, at 2:40 PM, Duncan Murdoch murdoch.dun...@gmail.com wrote: On 13-07-13 1:33 PM, Simon Zehnder wrote: Dear R-Users, I am writing a package using S4 classes. In the generic method plot I want to set the title for the plotting window as I will have several windows and window titles help the users to distinguish the graphics without putting a title into the plot itself (this can be done by users whenever they want) So I created a helper function .setDeviceTitle which I called after the plot has been done: .setDeviceTitle - function(title = title, dev = dev.cur()) { dev - names(dev)[1] ## check for OS ## if (dev == windows) { windows(title = title) } else if (dev == X11) { X11(title = title) } else { quartz(title = title) } } The result is a new device with the title in addition to the old. Is it possible to give a window a title after the plot has been done? If not: Before I plot the device I cannot know what device it will be, so I thought about a check via capabilities(): if (any(names(capabilities()) == X11)) { X11(title = title) } else if (any(names(capabilities)) == windows) { windows(title = title) } else { quartz(title = title) } I want to have a safe method, which works on each OS R can run. How would you solve the problem? Use dev.new() rather than picking a particular device. If all the possible devices support the title argument, then dev.new(title=title) will be fine. If you might need more customization (or want to protect against a user who chooses a device that doesn't have title as an argument), use getOption(device) to examine the call that will be used. Duncan Murdoch __ 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] Set window title for plot on any OS
Dear R-Users, I am writing a package using S4 classes. In the generic method plot I want to set the title for the plotting window as I will have several windows and window titles help the users to distinguish the graphics without putting a title into the plot itself (this can be done by users whenever they want) So I created a helper function .setDeviceTitle which I called after the plot has been done: .setDeviceTitle - function(title = title, dev = dev.cur()) { dev - names(dev)[1] ## check for OS ## if (dev == windows) { windows(title = title) } else if (dev == X11) { X11(title = title) } else { quartz(title = title) } } The result is a new device with the title in addition to the old. Is it possible to give a window a title after the plot has been done? If not: Before I plot the device I cannot know what device it will be, so I thought about a check via capabilities(): if (any(names(capabilities()) == X11)) { X11(title = title) } else if (any(names(capabilities)) == windows) { windows(title = title) } else { quartz(title = title) } I want to have a safe method, which works on each OS R can run. How would you solve the problem? Best Simon __ 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] Set window title for plot on any OS
Hi Duncan, thank you very much for your advice! That makes it all work. I check in addition for a title argument in the device via if (any(names(getOption(device)) == TRUE)) { dev.new(title = title) } Thanks again! Simon On Jul 13, 2013, at 2:40 PM, Duncan Murdoch murdoch.dun...@gmail.com wrote: On 13-07-13 1:33 PM, Simon Zehnder wrote: Dear R-Users, I am writing a package using S4 classes. In the generic method plot I want to set the title for the plotting window as I will have several windows and window titles help the users to distinguish the graphics without putting a title into the plot itself (this can be done by users whenever they want) So I created a helper function .setDeviceTitle which I called after the plot has been done: .setDeviceTitle - function(title = title, dev = dev.cur()) { dev - names(dev)[1] ## check for OS ## if (dev == windows) { windows(title = title) } else if (dev == X11) { X11(title = title) } else { quartz(title = title) } } The result is a new device with the title in addition to the old. Is it possible to give a window a title after the plot has been done? If not: Before I plot the device I cannot know what device it will be, so I thought about a check via capabilities(): if (any(names(capabilities()) == X11)) { X11(title = title) } else if (any(names(capabilities)) == windows) { windows(title = title) } else { quartz(title = title) } I want to have a safe method, which works on each OS R can run. How would you solve the problem? Use dev.new() rather than picking a particular device. If all the possible devices support the title argument, then dev.new(title=title) will be fine. If you might need more customization (or want to protect against a user who chooses a device that doesn't have title as an argument), use getOption(device) to examine the call that will be used. Duncan Murdoch __ 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] Optim seems not to work properly in foreach
Hi Ilai, after you sent this message I tried your code as well and it worked. As a result, I reconsidered the code written by me and of course also found the error in my function simulateEKOP. So for other people assuming errors or ill behavior in base functions: Forget it: these functions are tested this often, that the probabilitty is almost 1, the error is in your own code. Thanks for the help Simon On Jun 3, 2013, at 10:53 PM, ilai ke...@math.montana.edu wrote: On Mon, Jun 3, 2013 at 11:37 AM, Simon Zehnder szehn...@uni-bonn.de ... [Some not minimal, self contained, reproducible code]... Data simulation and thecreation of startpar works fine, but the parameters in res$par are always the start parameters. If I run the same commands directly on the shell I get in res$par the optimized parameters - only inside the foreach loop optim seems not to work. What could that be? Don't know, but but this makes me doubt it has anything to do with optim being inside foreach: fr - function(x) { x1 - x[1] ; x2 - x[2] 100 * (x2 - x1 * x1)^2 + (1 - x1)^2 } grr - function(x) { x1 - x[1] ; x2 - x[2] c(-400 * x1 * (x2 - x1 * x1) - 2 * (1 - x1) , 200 * (x2 - x1 * x1)) } library(doMC) registerDoMC(2) RNGkind(L'Ecuyer) set.seed(54321) foreach(i = 1:2) %do% { ret - foreach(j = 1:2) %do%{ strtpar - c(-2,2)+rnorm(2) optim(strtpar, fr, grr, method = L-BFGS-B,control=list(trace=TRUE))$par } ret } Also, wouldn't you want to register 4 cores by default if nesting 2 loops of 2 ? (to comment on the wisdom of doing so in terms of overhead is beyond my expertise) HTH Best Simon __ 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. [[alternative HTML version deleted]] __ 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] 3.0.1 update and compiler package
Hi Christoph, do you install from sources? Best Simon On Jun 4, 2013, at 10:41 AM, Christoph Knapp christoph.knap...@gmail.com wrote: Hi, reinstalling R did not help. It still will not update the same packages. I attached the terminal output. There were no errors while I was installing R again. I might have to remove all libraries and reinstall them all. Would you agree or do you think I should try something else first which is not as extreme. Thanks for your help Christoph On 03/06/13 20:32, Pascal Oettli wrote: My mistake, Regards, Pascal On 06/03/2013 04:37 PM, Uwe Ligges wrote: On 03.06.2013 07:19, Pascal Oettli wrote: Hi, How did you upgraded your version of R? From source or from a Linux package? Actually the new R installation is just broken. It simply has to be reinstalled carefully (watch for errors). Best, Uwe Ligges Regards, Pascal On 05/31/2013 11:33 AM, Christoph Knapp wrote: Hi, I recently updated to R 3.0.1. I'm running linux ubuntu 12.04. I realized that I have to update all the installed packages so I run update.packages(checkBuilt=TRUE) as described here http://www.r-bloggers.com/r-3-0-0-is-released-whats-new-and-how-to-upgrade/ . The first thing it did was telling me that it will not update several packages. When it was finished I used the warnings function to have a look at what did not work. See the list below. warnings() Warning messages: 1: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘boot’ had non-zero exit status 2: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘cluster’ had non-zero exit status 3: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘foreign’ had non-zero exit status 4: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘KernSmooth’ had non-zero exit status 5: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘MASS’ had non-zero exit status 6: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘Matrix’ had non-zero exit status 7: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘nlme’ had non-zero exit status 8: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘nnet’ had non-zero exit status 9: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘rpart’ had non-zero exit status 10: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘spatial’ had non-zero exit status 11: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘survival’ had non-zero exit status 12: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘class’ had non-zero exit status 13: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘epiR’ had non-zero exit status 14: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘gmodels’ had non-zero exit status 15: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘gplots’ had non-zero exit status 16: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘mgcv’ had non-zero exit status 17: In install.packages(update[instlib == l, Package], l, ... : installation of package ‘gregmisc’ had non-zero exit status I tried to reinstall them manually but this always failed because of a package dependency to the compiler package. Now, if I try to install the compiler package it tells me. install.packages(compiler) Installing package into ‘/home/christoph/R/x86_64-pc-linux-gnu-library/3.0’ (as ‘lib’ is unspecified) Warning message: package ‘compiler’ is not available (for R version 3.0.1) The last line also came up all the time when the packages were updated Doing a bit of research does not deliver much only that the compiler package was included into R at version 2.13.0 (http://dirk.eddelbuettel.com/blog/2011/04/12/). Most of those packages which do not work any more are pretty important for some of my scripts and I would not even know what packages replace the packages above. Would anyone know how to fix this? Regards __ 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-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] 3.0.1 update and compiler package
Hi Christoph, so what you could try to make sure it can run easily on your system is to install it from sources: First check if you have at least gcc version 4.7.2: gcc --version Open a terminal and type wget ftp://ftp.stat.math.ethz.ch/Software/R/R-patched.tar.gz tar -zxvf R-patched.tar.gz mkdir build cd build ../R-patched/configure --prefix=/usr/opt/R --with-blas='-framework libVec' --with-lapack Check the output at the end of the configuration process. If everything looks fine start the build: make -j4 This should run without errors. Then check if all works fine: make check If this is the case you can install it (up to here no installation have been done: so the folder /usr/opt/R does not exist yet): sudo make install Start R via the shell: /usr/opt/R/bin/R Install your packages and see what happens. If all is fine, add /usr/opt/R to your PATH variable, if not just delete the folder /usr/opt/R. It all worked perfectly on my scientific Linux system. Best Simon On Jun 4, 2013, at 11:50 AM, Christoph Knapp christoph.knap...@gmail.com wrote: Usually sudo apt-get install ... Sometimes synaptic package manager. I'm using the http://cran.ma.imperial.ac.uk/bin/linux/ubuntu/precise/ http://cran.stat.ucla.edu/bin/linux/ubuntu/precise repositories. Not sure what the differences are. I had some problems before when I installed packages from within Rstudio. It would install them in a different directory than apt-get would and it would tell me that those packages were not installed when I used the terminal to start R. Drove me crazy for a while because I could use them in Rstudio so I was sure they were installed. Do you think this could also cause this? Christoph On 04/06/13 21:13, Simon Zehnder wrote: Hi Christoph, do you install from sources? Best Simon On Jun 4, 2013, at 10:41 AM, Christoph Knapp christoph.knap...@gmail.com wrote: Hi, reinstalling R did not help. It still will not update the same packages. I attached the terminal output. There were no errors while I was installing R again. I might have to remove all libraries and reinstall them all. Would you agree or do you think I should try something else first which is not as extreme. Thanks for your help Christoph On 03/06/13 20:32, Pascal Oettli wrote: My mistake, Regards, Pascal On 06/03/2013 04:37 PM, Uwe Ligges wrote: On 03.06.2013 07:19, Pascal Oettli wrote: Hi, How did you upgraded your version of R? From source or from a Linux package? Actually the new R installation is just broken. It simply has to be reinstalled carefully (watch for errors). Best, Uwe Ligges Regards, Pascal On 05/31/2013 11:33 AM, Christoph Knapp wrote: Hi, I recently updated to R 3.0.1. I'm running linux ubuntu 12.04. I realized that I have to update all the installed packages so I run update.packages(checkBuilt=TRUE) as described here http://www.r-bloggers.com/r-3-0-0-is-released-whats-new-and-how-to-upgrade/ . The first thing it did was telling me that it will not update several packages. When it was finished I used the warnings function to have a look at what did not work. See the list below. warnings() Warning messages: 1: In install.packages(update[instlib == l, Package], l, ... : installation of package boot had non-zero exit status 2: In install.packages(update[instlib == l, Package], l, ... : installation of package cluster had non-zero exit status 3: In install.packages(update[instlib == l, Package], l, ... : installation of package foreign had non-zero exit status 4: In install.packages(update[instlib == l, Package], l, ... : installation of package KernSmooth had non-zero exit status 5: In install.packages(update[instlib == l, Package], l, ... : installation of package MASS had non-zero exit status 6: In install.packages(update[instlib == l, Package], l, ... : installation of package Matrix had non-zero exit status 7: In install.packages(update[instlib == l, Package], l, ... : installation of package nlme had non-zero exit status 8: In install.packages(update[instlib == l, Package], l, ... : installation of package nnet had non-zero exit status 9: In install.packages(update[instlib == l, Package], l, ... : installation of package rpart had non-zero exit status 10: In install.packages(update[instlib == l, Package], l, ... : installation of package spatial had non-zero exit status 11: In install.packages(update[instlib == l, Package], l, ... : installation of package survival had non-zero exit status 12: In install.packages(update[instlib == l, Package], l, ... : installation of package class had non-zero exit status 13: In install.packages(update[instlib == l, Package], l, ... : installation of package epiR had non-zero exit status 14: In install.packages(update[instlib == l, Package], l, ... : installation of package gmodels had non-zero exit status 15
[R] Optim seems not to work properly in foreach
Hi guys, I am working now for several years in R and I would say I manage things pretty easily, but with the foreach loop I have my problems. I call for a simulation a double foreach loop and this works fine. Inside the second loop (which I plan to parallelize later on) I call Rs optim-function: simulateML - function(sim.it = 250, input.path, output.path, T = 390, ncore = 2) { ## load packages ## library(mmstruct) library(doMC) library(foreach) ## read input parameters ## input - read.csv(input.path, header = FALSE) ## create container to store results ## results - data.frame(NA, nrow = NROW(input) * sim.it, ncol = 12) ## initialize the parallel backend ## registerDoMC(ncore) result.list - foreach(i = 1:2) %do% { input - as.matrix(input) input.row - input[i, ] list - foreach(i = 1:2, .combine = rbind, .packages = c(mmstruct, stats)) %do% { data - simulateEKOP(size = input.row[1], alpha = input.row[2], epsilon = input.row[1], delta = input.row[4], mu = input.row[5], T = T) data - data[,4] ## create start values ## tmp - mean(data)/T startpar - c(0, tmp * 0.75/2, tmp * 0.25/2) ## set options for optimization ## optim_fn- computeKokotLik optim_method- L-BFGS-B optim_lower - c(-1e+7, 0, 0) optim_upper - rep(1e+7, 3) optim_fnscale - -1 optim_maxit - 200 optim_ctrl - list(fnscale = optim_fnscale, maxit = optim_maxit) ## start optimization ## res - optim(par = startpar, fn = optim_fn, data = data, T = T, methodLik = approx, method = optim_method, lower = optim_lower, upper = optim_upper, control = optim_ctrl, hessian = TRUE) res$par } print(list) } } Data simulation and thecreation of startpar works fine, but the parameters in res$par are always the start parameters. If I run the same commands directly on the shell I get in res$par the optimized parameters - only inside the foreach loop optim seems not to work. What could that be? Best Simon __ 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] Count function calls
Dear Giovanni, apologize for this late reply! I was testing and reading a lot of stuff. I tried your suggestions and the problem of singularity in the regressor cross product vanishes when using the Group Mean function 'pgm' instead of 'pvcm'. Nevertheless, I found the collinearity in the regressors and 'pvcm' ran further, until I got another error. This time a little different but still directed towards the same area: 'Error in solve.default(crossprod(X[[i]])[!coefna[i, ], !coefna[i, ]])) : system is computationally singular: reciprocal condition number = 1.86131e-17' Following the description of R's function 'rcond', the reciprocal condition number measures how close a matrix is to be rank deficient. So in this case one regressor crossproduct seems to be sufficiently close to singularity, such that inversion becomes impossible. I attached you a simple subsample. If you let run the following commands on it: require(plm) data - read.csv(path, sep = ,, header = TRUE) pdata - plm.data(data, index = c(gvkey, fyear)) debug(plm:::pvcm) model.pvcm - pvcm(LOGDXSGA~LOGDSALE + DECRLOGDSALE + DECRLOGDSALEWDAVG, data = pdata, model = random) .. go until Browse[2] debug: ml - split(data cond) Then type: A - as.matrix(ml$'25062'[, 2:4]) XX - t(A) %*% A qr(XX)$rank [1] 2Aha, the rank is only 2! rcond(XX) [1] 9.339817e-16Alright, the matrix XX is pretty near to singularity, even it is not really singular! Let us have a look at it: A[, 2]/A[, 3] . A[, 2] is almost a multiple of A[, 3], if the latter is multiplied by -6.231979 or -6.231331 . at least it suffices to let the rank shrink to 2. The thing is, the same regression goes through in the Stata package. As I would say, I am an R fanatist, I would like to know, why it runs in Stata and not in R.and that one can let it run in R as well. Different number formats/precision? Maybe it suffices in such a case to enforce full rank of the crossproduct by enforcing positive definitess, for instance via the function 'nearPD' in the R package 'Matrix': Try library(Matrix) qr(nearPD(XX)$mat)$rank [1] 3 Yey!. :) Tell me your opinion Giovanni! Give me a little help here please! Best Simon P.S. Thank your for this fantastic package making panel data estimation possible in the R environment! On Feb 27, 2013, at 1:40 PM, Millo Giovanni giovanni_mi...@generali.com wrote: Hello. Another thing you may want to do depends on whether you are using model=within (the default) or model=random. In the first case, pvcm() estimates separate regressions, so you just need to loop lm() on individual indices to spot where it fails. In the second case, what you may want to do is try a similar estimator: pmg(..., model=mg), which is an unweighted version of Swamy's estimator in pvcm() and a simpler function to read and modify, possibly using the global assignment operator '-' as already suggested by William to output diagnostics. Best, Giovanni -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Simon Zehnder Sent: Tuesday, February 26, 2013 2:53 AM To: r-help@r-project.org help Subject: [R] Count function calls Dear R-users, I have the following problem: I am running the function 'pvcm' from the 'plm' Panel Data package. Inside this function 'solve' is called and gives for a certain individual data series an exception because of singularity. I would like to know which individual data series causes this error. I tried to debug it, but this is truly painful, as solve is called inside of 'lapply' and there are over 5,000 individual data series in the panel. Now, what I would like to do is to count the calls to 'solve' inside the function, so I can see, where the function throws the exception. I tried to use 'trace' with a count variable, but I have no clue how to define a global variable to be used by trace and updated at every call.is there another approach? Best Simon __ Ai sensi del D.Lgs. 196/2003 si precisa che le informazioni contenute in questo messaggio sono riservate ed a uso esclusivo del destinatario. Qualora il messaggio in parola Le fosse pervenuto per errore, La invitiamo ad eliminarlo senza copiarlo e a non inoltrarlo a terzi, dandocene gentilmente comunicazione. Grazie. Pursuant to Legislative Decree No. 196/2003, you are hereby informed that this message contains confidential information intended only for the use of the addressee. If you are not the addressee, and have received this message by mistake, please delete it and immediately notify us. You may not copy or disseminate this message to anyone. Thank you. __ R-help@r-project.org mailing list https://stat.ethz.ch
[R] S4-classes: Assignment of values to slots by reference
Dear R-users, I am working on a project that uses S4-classes. At some point I encountered the problem - well known to R - that I have to pass 3 different objects to a function, that should modify several slots of them and of course there is no passing by reference in R. Then I read this thread by Steve Lianoglou: https://stat.ethz.ch/pipermail/r-help/2010-August/250468.html, which offers from viewpoint of OOP an elegant solution. But I do not use an Object Method, but a regular function, pass in several objects and modify them. Here is an example with 1 Object using Steve's approach: setClass(example, representation( par = matrix, .cache = environment) ) setMethod(initialize, example, function(.Object, ..., .cache = new.env()) { callNextMethod(.Object, .cache = .cache, ...) } ) ex - new(example, par = matrix()) fmodify - function(O1) { pm - mean(rnorm(100)) pm - matrix(pm) O1@.cache$par - pm } fmodify(ex) So far so good. Everything worked nicely. But of course the value computed in the function 'fmodify' is now in the slot ex@.cache$par: ex@.cache$par and not in ex@par. Is there a possibility to get it into ex@par? Best Simon __ 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] Count function calls
Dear WIlliam, thank you very much for this very useful information! I learned something new! Enjoy your day! Simon On Feb 27, 2013, at 1:40 AM, William Dunlap wdun...@tibco.com wrote: This is where - is helpful: N - 0 ; trace(solve, quote(N - N + 1), print=FALSE) Tracing function solve in package base [1] solve lapply(3:0, function(i)solve(diag(i,3), 1:3)) Error in solve.default(diag(i, 3), 1:3) : Lapack routine dgesv: system is exactly singular: U[1,1] = 0 N [1] 4 You can also set options(error=recover) to look at the state of things when the error occurs. Bill Dunlap Spotfire, TIBCO Software wdunlap tibco.com -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Simon Zehnder Sent: Tuesday, February 26, 2013 2:53 AM To: r-help@r-project.org help Subject: [R] Count function calls Dear R-users, I have the following problem: I am running the function 'pvcm' from the 'plm' Panel Data package. Inside this function 'solve' is called and gives for a certain individual data series an exception because of singularity. I would like to know which individual data series causes this error. I tried to debug it, but this is truly painful, as solve is called inside of 'lapply' and there are over 5,000 individual data series in the panel. Now, what I would like to do is to count the calls to 'solve' inside the function, so I can see, where the function throws the exception. I tried to use 'trace' with a count variable, but I have no clue how to define a global variable to be used by trace and updated at every call.is there another approach? Best Simon __ 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-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] Count function calls
Dear R-users, I have the following problem: I am running the function 'pvcm' from the 'plm' Panel Data package. Inside this function 'solve' is called and gives for a certain individual data series an exception because of singularity. I would like to know which individual data series causes this error. I tried to debug it, but this is truly painful, as solve is called inside of 'lapply' and there are over 5,000 individual data series in the panel. Now, what I would like to do is to count the calls to 'solve' inside the function, so I can see, where the function throws the exception. I tried to use 'trace' with a count variable, but I have no clue how to define a global variable to be used by trace and updated at every call.is there another approach? Best Simon __ 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] Armadillo error in R extension
Dear Michael, during the last days of programming I came across RcppArmadillo (Rcpp) and for my kind of work this is indeed pretty interesting. So, I will take a closer look today what it is about and how it works. As a C++ programmer though, I am still interested why compilation of my package fails at this point of code. I hope I get a little more insight on the Rcpp-devel list suggested by Dirk! Thanks for your comment! Simon On Feb 2, 2013, at 12:14 AM, Michael Weylandt michael.weyla...@gmail.com wrote: Look up the terribly wonderful RcppArmadillo package. MW On Feb 1, 2013, at 8:38 PM, Simon Zehnder szehn...@uni-bonn.de wrote: Is there anyway with some experience in using armadillo in R C++ extensions? My problem is the following: I programmed a function in a header looking like #include armadillo inline arma::vec foo(input) { ... do something return an arma::vec object } compiling this via R CMD INSTALL packagename (PKG_CXXFLAGS = -I/folder/of/armadillo and armadillo_bits in my package) I get the following error error: expected initializer before 'foo' The exact line and word is where inline arma::vec ends. Does anyone know what this error could be? Best Simon __ 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-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] Armadillo error in R extension
Is there anyway with some experience in using armadillo in R C++ extensions? My problem is the following: I programmed a function in a header looking like #include armadillo inline arma::vec foo(input) { ... do something return an arma::vec object } compiling this via R CMD INSTALL packagename (PKG_CXXFLAGS = -I/folder/of/armadillo and armadillo_bits in my package) I get the following error error: expected initializer before 'foo' The exact line and word is where inline arma::vec ends. Does anyone know what this error could be? Best Simon __ 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] Modify objects in function
Dear R community, I do know, that an R function is constructing a copy of any object passed as argument into a function. I program on a larger S4 project for a package, and I arrived at a point where I have to think a little harder on implementation style (especially to spare users complex object handling). I have a function foo(), taking as input arguments two S4 objects of different class type foo - function(o1, o2) { o1@att1 - producesomething() o2@att2 - producesomethingelse() } Of course, this functions does not change the objects in the global environment. Now I have two choices 1. Change the objects and return a list with both objects: foo - function(o1, o2) { o1@att1 - producesomething() o2@att2 - producesomethingelse() l - list(O1 = o1, O2 = o2) return(l) } This is cumbersome for users, as they have then to assign the objects inside the returned list to the symbols used in the global environment. But it is the way intended by R. 2. Change the objects of the global environment inside the function: foo - function(o1, o2) { o1@att1 - producesomething() o2@att2 - producesomethingelse() assign(o1, o1, .GlobalEnv) assign(o2, o2, .GlobalEnv) } This is for users very practical, as the symbols used before refer now to objects with changed attributes and can be used immediately. BUT: It is not intended to be done like this. I spared any solution using one function for every single object type, as this becomes even more cumbersome than choice 1 above, in my opinion. What is your opinion on the trade-off between user-friendly handling of objects and R-intended programming style? Do I miss another choice, combining both? Best Simon __ 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] Modify objects in function
Dear Barry, thank you very much for this information. This looks pretty interesting! Best Simon On Jan 31, 2013, at 4:09 PM, Barry Rowlingson b.rowling...@lancaster.ac.uk wrote: it lets you do: (a~b~c) = foo() Mistook. should be: (a~b~c) %=% foo() because it defines the %=% operator. Barry __ 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] Modify objects in function
Hi Barry, this actually a good idea, to put them together! Probably even creating an object containing both of them. Haven't thought about it before. Best Simon On Jan 31, 2013, at 3:49 PM, Barry Rowlingson b.rowling...@lancaster.ac.uk wrote: On Thu, Jan 31, 2013 at 11:52 AM, Simon Zehnder szehn...@uni-bonn.de wrote: Dear R community, I do know, that an R function is constructing a copy of any object passed as argument into a function. I program on a larger S4 project for a package, and I arrived at a point where I have to think a little harder on implementation style (especially to spare users complex object handling). I have a function foo(), taking as input arguments two S4 objects of different class type foo - function(o1, o2) { o1@att1 - producesomething() o2@att2 - producesomethingelse() } Of course, this functions does not change the objects in the global environment. Now I have two choices 1. Change the objects and return a list with both objects: foo - function(o1, o2) { o1@att1 - producesomething() o2@att2 - producesomethingelse() l - list(O1 = o1, O2 = o2) return(l) } This is cumbersome for users, as they have then to assign the objects inside the returned list to the symbols used in the global environment. But it is the way intended by R. By cumbersome you mean the following (anti-?) pattern has to be used, for example in a loop: o1 = something() o2 = somethingelse() o12 = foo(o1,o2) o1 = o12$o1 o2 = o12$o2 so that at the end of it you have a modified o1 and o2? I suggest that if you have a function that modifies more than one of its arguments, then there is a strong case for making those arguments a single argument. So that you'd do: foo = function(o12){ o12$o1@att=bar() o12$o2@att=baz() return(o12) } o12 = list(something(), somethingelse()) o12 = foo(o12) and there you are, a modified o12 object. no packing/unpacking needed. I suspect that either your o1 and o2 are so closely related that you should pack them in a list (or ideally, an object) or differently related such that you should treat them separately and not tweak both of them within the same function. That approach is binding behaviour very tightly to two structures, and probably won't give you very debuggable modular code. B __ 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] Irreproducible exception in R extension
Dear R-Users, I am having some trouble running an R extension on our cluster (linux). I call C++ code in which I use a) the Scythe Statistical Library and b) OpenMP. Most of the jobs run without a problem, but some arbitrary jobs throw an exception of the kind printed below while running in a parallel loop. The behavior is to me not reproducible, although at every run of 729 Jobs it at least happens once. I already tried to enforce the error by using the Intel Compiler instead of the GCC 4.7.2 and using the Options '-checkinit' for possibly uninitialized values and '-ftrapuv' for initialization of stack variables to unusual values to aid error detection. But I cannot enforce the exception, which makes it impossible to debug. I already googled this error, but there only 2-3 small forum discussions about it, and so I hope I find some suggestions here. 1. What do you think this error comes from? 2. What do you advice me to do to detect the error? I am excited for your answers. Simon *** glibc detected *** /home/simo_109//opt/R/lib64/R/bin/exec/R: double free or corruption (!prev): 0x01355610 *** === Backtrace: = /lib64/libc.so.6[0x3845a75916] /lib64/libc.so.6[0x3845a78443] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_ZN6scythe9DataBlockIdE10deallocateEv+0x7b)[0x7f4645d5d6db] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_ZN6scythe9DataBlockIdED1Ev+0x47)[0x7f4645d5d5bb] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_ZN6scythe18DataBlockReferenceIdE17withdrawReferenceEv+0xaa)[0x7f4645d83442] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_ZN6scythe18DataBlockReferenceIdED1Ev+0x94)[0x7f4645d82cc8] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_ZN6scythe18DataBlockReferenceIdED2Ev+0x47)[0x7f4645d82d9f] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_ZN6scythe6MatrixIdLNS_12matrix_orderE0ELNS_12matrix_styleE1EED1Ev+0x5d)[0x7f4645d8b26d] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_Z11MCgmmS_implIN6scythe8mersenneEEvRNS0_3rngIT_EERP7SEXPRECS8_jjRNS0_6MatrixIdLNS0_12matrix_orderE0ELNS0_12matrix_styleE0EEERKSC_jS8_S8_+0x478f)[0x7f4645d6d6af] /opt/intel/Compiler/12.1/6.361/rwthlnk/compiler/lib/intel64/libiomp5.so(__kmp_invoke_microtask+0x93)[0x7f4645aabee3] === Memory map: 0040-00401000 r-xp 00:21 32550383 /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/bin/exec/R 0060-00601000 rw-p 00:21 32550383 /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/bin/exec/R 00601000-029cc000 rw-p 00:00 0 [heap] 369ac0-369ac3a000 r-xp 08:01 263771 /lib64/libreadline.so.6.0 369ac3a000-369ae3a000 ---p 0003a000 08:01 263771 /lib64/libreadline.so.6.0 369ae3a000-369ae42000 rw-p 0003a000 08:01 263771 /lib64/libreadline.so.6.0 369ae42000-369ae43000 rw-p 00:00 0 369bc0-369bc16000 r-xp 08:01 263757 /lib64/libgcc_s-4.4.6-20120305.so.1 369bc16000-369be15000 ---p 00016000 08:01 263757 /lib64/libgcc_s-4.4.6-20120305.so.1 369be15000-369be16000 rw-p 00015000 08:01 263757 /lib64/libgcc_s-4.4.6-20120305.so.1 369c00-369c0e8000 r-xp 08:01 656025 /usr/lib64/libstdc++.so.6.0.13 369c0e8000-369c2e8000 ---p 000e8000 08:01 656025 /usr/lib64/libstdc++.so.6.0.13 369c2e8000-369c2ef000 r--p 000e8000 08:01 656025 /usr/lib64/libstdc++.so.6.0.13 369c2ef000-369c2f1000 rw-p 000ef000 08:01 656025 /usr/lib64/libstdc++.so.6.0.13 369c2f1000-369c306000 rw-p 00:00 0 369d00-369d00d000 r-xp 08:01 671072 /usr/lib64/libgomp.so.1.0.0 369d00d000-369d20c000 ---p d000 08:01 671072 /usr/lib64/libgomp.so.1.0.0 369d20c000-369d20d000 rw-p c000 08:01 671072 /usr/lib64/libgomp.so.1.0.0 369fa0-369faf r-xp 08:01 660772 /usr/lib64/libgfortran.so.3.0.0 369faf-369fcef000 ---p 000f 08:01 660772 /usr/lib64/libgfortran.so.3.0.0 369fcef000-369fcf1000 rw-p 000ef000 08:01 660772 /usr/lib64/libgfortran.so.3.0.0 369fcf1000-369fcf2000 rw-p 00:00 0 36a360-36a373f000 r-xp 08:01 675210 /usr/lib64/libicuuc.so.42.1 36a373f000-36a393f000 ---p 0013f000 08:01 675210 /usr/lib64/libicuuc.so.42.1 36a393f000-36a395 rw-p 0013f000 08:01 675210 /usr/lib64/libicuuc.so.42.1 36a395-36a3952000 rw-p 00:00 0 36a3a0-36a3b88000 r-xp
[R] Irreproducible exception in R extension
Dear R-Users, I am having some trouble running an R extension on our cluster (linux). I call C++ code in which I use a) the Scythe Statistical Library and b) OpenMP. Most of the jobs run without a problem, but some arbitrary jobs throw an exception of the kind printed below while running in a parallel loop. The behavior is to me not reproducible, although at every run of 729 Jobs it at least happens once. I already tried to enforce the error by using the Intel Compiler instead of the GCC 4.7.2 and using the Options '-checkinit' for possibly uninitialized values and '-ftrapuv' for initialization of stack variables to unusual values to aid error detection. But I cannot enforce the exception, which makes it impossible to debug. I already googled this error, but there only 2-3 small forum discussions about it, and so I hope I find some suggestions here. 1. What do you think this error comes from? 2. What do you advice me to do to detect the error? I am excited for your answers. Simon *** glibc detected *** /home/simo_109//opt/R/lib64/R/bin/exec/R: double free or corruption (!prev): 0x01355610 *** === Backtrace: = /lib64/libc.so.6[0x3845a75916] /lib64/libc.so.6[0x3845a78443] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_ZN6scythe9DataBlockIdE10deallocateEv+0x7b)[0x7f4645d5d6db] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_ZN6scythe9DataBlockIdED1Ev+0x47)[0x7f4645d5d5bb] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_ZN6scythe18DataBlockReferenceIdE17withdrawReferenceEv+0xaa)[0x7f4645d83442] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_ZN6scythe18DataBlockReferenceIdED1Ev+0x94)[0x7f4645d82cc8] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_ZN6scythe18DataBlockReferenceIdED2Ev+0x47)[0x7f4645d82d9f] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_ZN6scythe6MatrixIdLNS_12matrix_orderE0ELNS_12matrix_styleE1EED1Ev+0x5d)[0x7f4645d8b26d] /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/library/MCpkg/libs/MCpkg.so(_Z11MCgmmS_implIN6scythe8mersenneEEvRNS0_3rngIT_EERP7SEXPRECS8_jjRNS0_6MatrixIdLNS0_12matrix_orderE0ELNS0_12matrix_styleE0EEERKSC_jS8_S8_+0x478f)[0x7f4645d6d6af] /opt/intel/Compiler/12.1/6.361/rwthlnk/compiler/lib/intel64/libiomp5.so(__kmp_invoke_microtask+0x93)[0x7f4645aabee3] === Memory map: 0040-00401000 r-xp 00:21 32550383 /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/bin/exec/R 0060-00601000 rw-p 00:21 32550383 /rwthfs/rz/cluster/home/simo_109/opt/R/lib64/R/bin/exec/R 00601000-029cc000 rw-p 00:00 0 [heap] 369ac0-369ac3a000 r-xp 08:01 263771 /lib64/libreadline.so.6.0 369ac3a000-369ae3a000 ---p 0003a000 08:01 263771 /lib64/libreadline.so.6.0 369ae3a000-369ae42000 rw-p 0003a000 08:01 263771 /lib64/libreadline.so.6.0 369ae42000-369ae43000 rw-p 00:00 0 369bc0-369bc16000 r-xp 08:01 263757 /lib64/libgcc_s-4.4.6-20120305.so.1 369bc16000-369be15000 ---p 00016000 08:01 263757 /lib64/libgcc_s-4.4.6-20120305.so.1 369be15000-369be16000 rw-p 00015000 08:01 263757 /lib64/libgcc_s-4.4.6-20120305.so.1 369c00-369c0e8000 r-xp 08:01 656025 /usr/lib64/libstdc++.so.6.0.13 369c0e8000-369c2e8000 ---p 000e8000 08:01 656025 /usr/lib64/libstdc++.so.6.0.13 369c2e8000-369c2ef000 r--p 000e8000 08:01 656025 /usr/lib64/libstdc++.so.6.0.13 369c2ef000-369c2f1000 rw-p 000ef000 08:01 656025 /usr/lib64/libstdc++.so.6.0.13 369c2f1000-369c306000 rw-p 00:00 0 369d00-369d00d000 r-xp 08:01 671072 /usr/lib64/libgomp.so.1.0.0 369d00d000-369d20c000 ---p d000 08:01 671072 /usr/lib64/libgomp.so.1.0.0 369d20c000-369d20d000 rw-p c000 08:01 671072 /usr/lib64/libgomp.so.1.0.0 369fa0-369faf r-xp 08:01 660772 /usr/lib64/libgfortran.so.3.0.0 369faf-369fcef000 ---p 000f 08:01 660772 /usr/lib64/libgfortran.so.3.0.0 369fcef000-369fcf1000 rw-p 000ef000 08:01 660772 /usr/lib64/libgfortran.so.3.0.0 369fcf1000-369fcf2000 rw-p 00:00 0 36a360-36a373f000 r-xp 08:01 675210 /usr/lib64/libicuuc.so.42.1 36a373f000-36a393f000 ---p 0013f000 08:01 675210 /usr/lib64/libicuuc.so.42.1 36a393f000-36a395 rw-p 0013f000 08:01 675210 /usr/lib64/libicuuc.so.42.1 36a395-36a3952000 rw-p 00:00 0 36a3a0-36a3b88000 r-xp
Re: [R] Generating an autocorrelated binary variable
Thank you very much for your answer Rolf. It helped. I try to simulate a trade indicator model from market microstructure, where the 1 or -1 indicate a buyer or seller initiated trade respectively. I use a Gaussian copula for simulation, so I can put in some correlation if I want to. So I generate my multivariate Gaussian sample and transform it to a uniform sample so I can use whatever marginal distribution I want to. Then I generate the autocorrelated binary by beginning with a transformation on the first observation and then by using the uniform sample to compare it to 0.6. Afterwards I have a 0.6 autocorrelated binary variable: sampleCop - function(n = 1000, rho = 0.6) { require(splus2R) mvrs - rmvnorm(n + 1, mean = rep(0, 3), cov = diag(3)) pmvrs - pnorm(mvrs, 0, 1) var1 - matrix(0, nrow = n + 1, ncol = 1) var1[1] - qbinom(pmvrs[1, 1], 1, 0.5) if(var1[1] == 0) var1[1] - -1 for(i in 1:(nrow(pmvrs) - 1)) { if(pmvrs[i + 1, 1] = rho) var1[i + 1] - var1[i] else var1[i + 1] - var1[i] * (-1) } sample - matrix(0, nrow = n, ncol = 4) sample[, 2] - var1[1:nrow(var1) - 1] sample[, 3] - var1[2:nrow(var1)] sample[, 4] - qnorm(pmvrs[1:nrow(var1) - 1, 2], 0, 1, 1, 0) + qnorm(pmvrs[1:nrow(var1) - 1, 3], 0, 1, 1, 0) - qnorm(pmvrs[2:nrow(var1), 3], 0, 1, 1, 0) sample[, 1] - sample[,2] * (0.015 + 0.03) - sample[,3] * (0.015 + 0.2*0.03) + sample[,4] sample } It is interesting though, that if I run my GMM-model on it it needs very high numbers of observations to converge to the real parameters. Using the gmm-function gmmfun - function(theta, data) { res = data[,1] - (theta[1] + theta[2]) * data[,2] + (theta[1] + theta[3] * theta[2]) * data[,3] moments = matrix(0, nrow=nrow(data), ncol=3) moments[,1] = data[,2] * data[,3] - data[,2]^2 * theta[3] moments[,2] = res * data[,2] moments[,3] = res * data[,3] return(moments) } And the objective function obj_fun - function(theta, data) { moments - gmmfun(theta, data) moments_sum - colSums(moments) val - moments_sum %*% diag(3) %*% moments_sum * 1.0 / nrow(data) val } Running then the command optim(fn = obj_fun, par = c(0.05, 0.05, .1), method = BFGS) gives a result with the first parameter fairly high away from its true value (true values for the three parameters were (0.015, 0.03, 0.2)). I wonder if this is due to the autocorrelation in the binary variable (btw. I also tried here lm, as the third parameter can be estimated via the correlation, but it remains highly dependent on the number of observations). Any suggestions here? Best regards Simon On Sep 28, 2012, at 5:02 AM, Rolf Turner rolf.tur...@xtra.co.nz wrote: I have no idea what your code is doing, nor why you want correlated binary variables. Correlation makes little or no sense in the context of binary random variables --- or more generally in the context of discrete random variables. Be that as it may, it is an easy calculation to show that if X and Y are binary random variables both with success probability of 0.5 then cor(X,Y) = 0.2 if and only if Pr(X=1 | Y = 1) = 0.6. So just generate X and Y using that fact: set.seed(42) X - numeric(1000) Y - numeric(1000) for(i in 1:1000) { Y[i] - rbinom(1,1,0.5) X[i] - if(Y[i]==1) rbinom(1,1,0.6) else rbinom(1,1,0.4) } # Check: cor(X,Y) # Get 0.2012336 Looks about right. Note that the sample proportions are 0.484 and 0.485 for X and Y respectively. These values do not differ significantly from 0.5. cheers, Rolf Turner On 28/09/12 08:26, Simon Zehnder wrote: Hi R-fellows, I am trying to simulate a multivariate correlated sample via the Gaussian copula method. One variable is a binary variable, that should be autocorrelated. The autocorrelation should be rho = 0.2. Furthermore, the overall probability to get either outcome of the binary variable should be 0.5. Below you can see the R code (I use for simplicity a diagonal matrix in rmvnorm even if it produces no correlated sample): sampleCop - function(n = 1000, rho = 0.2) { require(splus2R) mvrs - rmvnorm(n + 1, mean = rep(0, 3), cov = diag(3)) pmvrs - pnorm(mvrs, 0, 1) var1 - matrix(0, nrow = n + 1, ncol = 1) var1[1] - qbinom(pmvrs[1, 1], 1, 0.5) if(var1[1] == 0) var1[nrow(mvrs)] - -1 for(i in 1:(nrow(pmvrs) - 1)) { if(pmvrs[i + 1, 1] = rho) var1[i + 1] - var1[i] else var1[i + 1] - var1[i] * (-1) } sample - matrix(0, nrow = n, ncol = 4) sample[, 1] - var1[1:nrow(var1) - 1] sample[, 2] - var1[2:nrow(var1)] sample[, 3] - qnorm(pmvrs[1:nrow(var1) - 1, 2], 0, 1, 1, 0) sample[, 4] - qnorm(pmvrs[1:nrow(var1) - 1, 3], 0, 1, 1, 0) sample
[R] Autocorrelated binary variable
Hi R-fellows, I am trying to simulate a multivariate correlated sample via the Gaussian copula method. One variable is a binary variable, that should be autocorrelated. The autocorrelation should be rho = 0.2. Furthermore, the overall probability to get either outcome of the binary variable should be 0.5. Below you can see the R code (I use for simplicity a diagonal matrix in rmvnorm even if it produces no correlated sample): sampleCop - function(n = 1000, rho = 0.2) { require(splus2R) mvrs - rmvnorm(n + 1, mean = rep(0, 3), cov = diag(3)) pmvrs - pnorm(mvrs, 0, 1) var1 - matrix(0, nrow = n + 1, ncol = 1) var1[1] - qbinom(pmvrs[1, 1], 1, 0.5) if(var1[1] == 0) var1[nrow(mvrs)] - -1 for(i in 1:(nrow(pmvrs) - 1)) { if(pmvrs[i + 1, 1] = rho) var1[i + 1] - var1[i] else var1[i + 1] - var1[i] * (-1) } sample - matrix(0, nrow = n, ncol = 4) sample[, 1] - var1[1:nrow(var1) - 1] sample[, 2] - var1[2:nrow(var1)] sample[, 3] - qnorm(pmvrs[1:nrow(var1) - 1, 2], 0, 1, 1, 0) sample[, 4] - qnorm(pmvrs[1:nrow(var1) - 1, 3], 0, 1, 1, 0) sample } Now, the code is fine, everything compiles. But when I compute the autocorrelation of the binary variable, it is not 0.2, but 0.6. Does anyone know why this happens? Best Regards Simon [[alternative HTML version deleted]] __ 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] Generating an autocorrelated binary variable
Hi R-fellows, I am trying to simulate a multivariate correlated sample via the Gaussian copula method. One variable is a binary variable, that should be autocorrelated. The autocorrelation should be rho = 0.2. Furthermore, the overall probability to get either outcome of the binary variable should be 0.5. Below you can see the R code (I use for simplicity a diagonal matrix in rmvnorm even if it produces no correlated sample): sampleCop - function(n = 1000, rho = 0.2) { require(splus2R) mvrs - rmvnorm(n + 1, mean = rep(0, 3), cov = diag(3)) pmvrs - pnorm(mvrs, 0, 1) var1 - matrix(0, nrow = n + 1, ncol = 1) var1[1] - qbinom(pmvrs[1, 1], 1, 0.5) if(var1[1] == 0) var1[nrow(mvrs)] - -1 for(i in 1:(nrow(pmvrs) - 1)) { if(pmvrs[i + 1, 1] = rho) var1[i + 1] - var1[i] else var1[i + 1] - var1[i] * (-1) } sample - matrix(0, nrow = n, ncol = 4) sample[, 1] - var1[1:nrow(var1) - 1] sample[, 2] - var1[2:nrow(var1)] sample[, 3] - qnorm(pmvrs[1:nrow(var1) - 1, 2], 0, 1, 1, 0) sample[, 4] - qnorm(pmvrs[1:nrow(var1) - 1, 3], 0, 1, 1, 0) sample } Now, the code is fine, everything compiles. But when I compute the autocorrelation of the binary variable, it is not 0.2, but 0.6. Does anyone know why this happens? Best Regards Simon [[alternative HTML version deleted]] __ 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] openmp on R?
Dear R-Users, in the configuring file of R on the supercomputer I am using I found a line --disable-openmp (version R 2.13.1) Is there a possibility to enable this when I let run a BATCH file via R CMD? Thanks for comments! Simon __ 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] Problem with memory consuming algorithm
Hi guyz, I have serious problems with an algorithm I let run on a supercomputer: You find the functions under the following URLs: simuFunctionCaller: http://pastebin.com/6gw2fJFb calls Function simuFunctionNBM (http://pastebin.com/QeJDUnqx) after reading a csv-file ordered like the following: ,name,theta,std_error_theta,theta_phi,std_error_theta_phi,phi,N,sigma_epsilon,sigma_eta 31,25%,0.003044,2.28E-09,0.00742575,4.71E-09,0.00452525,68164.75,0.008706072,0.00774016 32,50%,0.004536,3.22E-09,0.01045,7.94E-09,0.0052762,110150,0.012387386,0.011169953 33,75%,0.007356625,1.24E-08,0.014249075,3.16E-08,0.007518,178840.5,0.018471461,0.016058424 Reading and running without an exception is no problem. But the algorithm needs hours to run and seems to use more memory than available: I have given a time limit of 48 h and a memory limit of 128 GB RAM on the supercomputer (it is a cluster, a BATCH System)…..nevertheless the computer shuts down the algorithm after a certain time because of memory problems. Does anyone of u guys see a problem in the algorithm? Especially with so much RAM and time available? I am very thankful for your suggestions! Simon __ 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] Exception in NeweyWest - Pre-Whitening necessary?
Hi guyz, I have run my algorithm in R (see http://pastebin.com/q84Tujfg) and got the following error: Error in ar.ols(x, aic = aic, order.max = order.max, na.action = na.action, : 'order.max' must be 'n.used' I am pretty sure, that the error comes from the NeweyWest function in line 45, as the NeweyWest function uses the ar.ols() function for pre whitening. Does anyone has an idea how to circumvent this? Can I just shut off pre whitening? I am thankful for every suggestion. Simon [[alternative HTML version deleted]] __ 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] Exception while using NeweyWest function with doMC
Hi Jay, first: thank u very much for your comments! U made some very important points clear. I tried immediately to write directly the sample function from trade-as.big.matrix(matrix(sample(c(1,-1), (N+1)*K, replace=TRUE),ncol=K), backingpath=backingpath, backingfile=trade.bin,descriptorfile=trade.desc) into the big matrix: trade-big.matrix(sample(c(1,-1), (10+1), replace=TRUE),nrow=(10+1), ncol=10, type=double,backingpath=/Users/simon/Documents/R/BigMTest/, backingfile=terminaltest.bin, descriptorfile=terminaltest.desc) But I either get only 1s: trade[,] [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [1,]111111111 1 [2,]111111111 1 [3,]111111111 1 [4,]111111111 1 [5,]111111111 1 [6,]111111111 1 [7,]111111111 1 [8,]111111111 1 [9,]111111111 1 [10,]111111111 1 [11,]111111111 1 or only -1s: trade[,] [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [1,] -1 -1 -1 -1 -1 -1 -1 -1 -1-1 [2,] -1 -1 -1 -1 -1 -1 -1 -1 -1-1 [3,] -1 -1 -1 -1 -1 -1 -1 -1 -1-1 [4,] -1 -1 -1 -1 -1 -1 -1 -1 -1-1 [5,] -1 -1 -1 -1 -1 -1 -1 -1 -1-1 [6,] -1 -1 -1 -1 -1 -1 -1 -1 -1-1 [7,] -1 -1 -1 -1 -1 -1 -1 -1 -1-1 [8,] -1 -1 -1 -1 -1 -1 -1 -1 -1-1 [9,] -1 -1 -1 -1 -1 -1 -1 -1 -1-1 [10,] -1 -1 -1 -1 -1 -1 -1 -1 -1-1 [11,] -1 -1 -1 -1 -1 -1 -1 -1 -1-1 Is there another possibility? In addition I found under ?as.big.matrix, for the second example the usage of matrix() inside of as.big.matrix(). But if I understood u correctly: usage is possible but does not save memory? I used the big.memory package because I always got - using simply matrix() - exceptions, telling me, that memory has reached its limits. After using big.memory all worked fine, BUT running http://pastebin.com/UxSkzrae and http://pastebin.com/MErGQsQd , I got this: http://pastebin.com/KrEncrSz. It seems, as there is a problem with memory allocation inside the underlying C-code, maybe a result from my matrix generation inside of the big matrix? Any suggestions? On Aug 29, 2011, at 6:24 PM, Jay Emerson wrote: Simon, Though we're please to see another use of bigmemory, it really isn't clear that it is gaining you anything in your example; anything like as.big.matrix(matrix(...)) still consumes full RAM for both the inner matrix() and the new big.matrix -- is the filebacking really necessary. It also doesn't appear that you are making use of shared memory, so I'm unsure what the gains are. However, I don't have any particular insight as to the subsequent problem with NeweyWest (which doesn't seem to be using the big.matrix objects). Jay -- Message: 32 Date: Sat, 27 Aug 2011 21:37:55 +0200 From: Simon Zehnder simon.zehn...@googlemail.com To: r-help@r-project.org Subject: [R] Exception while using NeweyWest function with doMC Message-ID: cagqvrp_gk+t0owbv1ste-y0zafmi9s_zwqrxyxugsui18ms...@mail.gmail.com Content-Type: text/plain Dear R users, I am using R right now for a simulation of a model that needs a lot of memory. Therefore I use the *bigmemory* package and - to make it faster - the *doMC* package. See my code posted on http://pastebin.com/dFRGdNrG snip - -- John W. Emerson (Jay) Associate Professor of Statistics Department of Statistics Yale University http://www.stat.yale.edu/~jay [[alternative HTML version deleted]] __ 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] Descriptive Stats from Data Frame
Hi Rich, I do not know what u really want, because it seems to me, u want to calculate the mean of all rows, where the chemical is Arsenic?? But try this to get a little more inside: mean(chemdata$quant[chemdata$param==Arsenic]) The vector chemdata[chemdata$param==Arsenic,] is a logical vector, returning TRUE for every row in which the variable param takes the value Arsenic. Try it in your R editor to see it and understand the R concept! If u now want to get all values of a certain column, given all values have Arsenic as param, u just write: chemdata$COLUMNNAME[chemdata$param==Arsenic] I do not know if this helps, as it seems to me, that Arsenic only occurs once in your frame….. Good luck Simon On Aug 30, 2011, at 11:00 PM, Rich Shepard wrote: I don't find how to do what I need to do in Dalgaard or 'R Cookbook', so I'm asking here. I have a data frame with water chemistry data and I want to start exploring these data. There are three factors (site, date, chemical) associated with each measurement. The data frame looks like this: summary(chemdata) site_id.sample_date.param.quant BC-0.5|1996-04-19|Arsenic|0.01 :1 BC-0.5|1996-04-19|Calcium|76.56 :1 BC-0.5|1996-04-19|Chloride|12 :1 BC-0.5|1996-04-19|Magnesium|43.23 :1 BC-0.5|1996-04-19|Sulfate|175 :1 BC-0.5|1996-04-19|Total Dissolved Solids|460:1 (Other) :14880 I want first to calculate (and plot) descriptive stats by chemical, ignoring site and date and telling R to ignore missing data. (Incorporating those factors will occur later.) What I have not been able to figure out is how to specify the command to, for example, calculate mean and sd for Arsenic. My floundering and thrashing includes attempts like these: mean(chemdata.param=Arsenic) Error in is.numeric(x) : 'x' is missing mean(chemdata.quant, param=Arsenic) Error in mean(chemdata.quant, param = Arsenic) : object 'chemdata.quant' not found mean(chemdata$quant, param=Arsenic) [1] NA Warning message: In mean.default(chemdata$quant, param = Arsenic) : argument is not numeric or logical: returning NA As a newcomer to R I've done a lot of reading, yet all the examples use nicely structured data to illustrate the point being made. I need to work with my data and learn how to specify columns and write correct commands for the analyses I need. Please point me in the right direction. Rich __ 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-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] Exception while using NeweyWest function with doMC
Dear R users, I am using R right now for a simulation of a model that needs a lot of memory. Therefore I use the *bigmemory* package and - to make it faster - the *doMC* package. See my code posted on http://pastebin.com/dFRGdNrG Now, if I use the foreach loop with the addon %do% (for sequential run) I have no problems at all - only here and there some singularities in regressor matrices which should be ok. BUT if I run the loop on multiple cores I get very often a bad exception. I have posted the exception on http://pastebin.com/eMWF4cu0 The exception comes from the NeweyWest function loaded within the sandwich library. I have no clue, what it want to say me and why it is so weirdly printed to the terminal. I am used to receive here and there errorsbut the messages never look like this. Does anyone have a useful answer for me, where to look for the cause of this weird error? Here some additional information: Hardware: MacBook Pro 2.66 GHz Intel Core Duo, 4 GB Memory 1067 MHz DDR3 Software System: Mac Os X Lion 10.7.1 (11B26) Software App: R64 version 2.11.1 run via Mac terminal I hope someone has a good suggestion! Thank u all! Simon [[alternative HTML version deleted]] __ 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.