Re: [R] How to change the value of a class slot
Ross Boylan wrote:
On Tue, 2005-06-07 at 08:36 +0200, Uwe Ligges wrote:
Ross Boylan wrote:
On Sat, Jun 04, 2005 at 07:08:56PM +0200, Martin Maechler wrote:
Ross nextPath - function(pm){ #pm is a CompletePathMaker
Ross[EMAIL PROTECTED] - [EMAIL PROTECTED](1)
Ross [etc]
If your nextPath function has 'pm' as its last statement it
will return the updated object, and if you call it
as
mypm - nextPath(mypm)
you are
1) updating mypm
2) in a proper S way (i.e. no cheating).
Regards,
Martin
Wow. This is almost the exact inverse of the usual object behavior,
in which only the class itself can update the slots (aka instance
variables). None of the methods of the class can update instances of
the class persistently without the help of outsiders, and only
outsiders can change the slot values.
(Yes, I realize that using the idiom you suggest of returning a new
object one can have only class methods actually fiddling with the
slots.)
The inability of a class method to change a class object without
outside help seems unfortunate.
It looks as if instances of class objects are best thought of as
immutable once created.
Obviously, there are many definition of object oriented programming,
and yours seems to be different from the S4 definition.
Yes. And though there are many definitions of object oriented (at
least, many implementations), I'd say the minimum requirement to be
object oriented is to have objects that encapsulate both state (instance
variables/slots) and behavior (methods).
S4 objects do not fully encapsulate state because they require outside
assistance to alter the state of the object (with the exception of
assignment operators). The smalltalker in me also gets nervous that
code outside the class can access the slots, but there are many object
systems that act that way.
The way in which names of methods of unrelated classes interfere with
each other seems a break-down of the encapsulation of behavior, though
the problem strictly is not with the behavior but just with the name.
To return to the concrete problem that got me started, if class
Specification defines a method likelihood taking as arguments instances
of class Specification, Path and Parameters, then it is awkward to
define a method likelihood for the class Model when that method has
arguments of class Model, Specification, data.frame, and vector,
particularly if different names for the formal arguments are desired.
(I think technically it could be done, but only in a very ugly
way--i.e., better to use different method names for the two classes).
I was going to answer your first question at first, but you have not
given enough details - in particular it was not clear to me why your
approach did not work.
I assumed that you are assigning the new object
again, which is the S way.
I wasn't, which is why it didn't work. I wanted the function to return
some other value than the object it was operating on.
You have to think about scoping rules and it
will be clear that the approach you are expecting is not a clean one in S.
Could you say a bit more about that?
I meant the following simple example (not related to any object oriented
programming from the S point of view, but maybe well from your point of
view?):
Say you want a function foo() that simply incements the argument:
a - 1
foo(a)
a # now is 2
But what happens if there is (more than) one a (in different
environments), but no a in the environment foo(a) is called from.
Which a do you want to change in this case? Seems to be rather dangerous.
Uwe Ligges
I had thought of the issue more in
terms of function calls in S being call by value, preventing updates to
the original arguments. So the issue isn't so much the scope of the
names of function arguments (that scope being limited to the function
body), but the properties of the thing they refer to (conceptually, a
copy of the argument, not the original).
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] matrix
hi, is it possible to create a matrix with one row and according to add a row? in fact, at present, I'm doing an algorithm which fill a matrix. On the web site of CRAN, the package basic is impossible to be load! can you tell me where i can found it! Thanks Sabine - [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] matrix
On Wed, 8 Jun 2005 10:37:15 +0200 (CEST) Navarre Sabine wrote: hi, is it possible to create a matrix with one row and according to add a row? in fact, at present, I'm doing an algorithm which fill a matrix. I think you are looking for ?rbind and might also be interested in looking at An Introduction to R or other introductory manuals. On the web site of CRAN, the package basic is impossible to be load! can you tell me where i can found it! You already asked this question yesterday (and received no reply because it doesn't make sense). There is no package with the name basic. There is only the base package which is part of every distribution of R (which I hope you have installed) and loaded at startup of R. Z __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] matrix
Achim Zeileis wrote: On Wed, 8 Jun 2005 10:37:15 +0200 (CEST) Navarre Sabine wrote: hi, is it possible to create a matrix with one row and according to add a row? in fact, at present, I'm doing an algorithm which fill a matrix. I think you are looking for ?rbind and might also be interested in looking at An Introduction to R or other introductory manuals. Achim is completely correct, I'd like to add: If you do know the size of the matrix before your algorithm starts iterating, you should not use rbind() but generate a fullsize matrix at first and replace row by row afterwards. Uwe Ligges On the web site of CRAN, the package basic is impossible to be load! can you tell me where i can found it! You already asked this question yesterday (and received no reply because it doesn't make sense). There is no package with the name basic. There is only the base package which is part of every distribution of R (which I hope you have installed) and loaded at startup of R. Z __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] matrix
Le 08.06.2005 10:37, Navarre Sabine a écrit : hi, is it possible to create a matrix with one row and according to add a row? in fact, at present, I'm doing an algorithm which fill a matrix. ?rbind When you know the final size of your matrix, it is better to create a full matrix and then fill it. On the web site of CRAN, the package basic is impossible to be load! can you tell me where i can found it! I don't see any package called 'basic' on CRAN. What is that package ? Thanks Sabine -- visit the R Graph Gallery : http://addictedtor.free.fr/graphiques ~ ~~ Romain FRANCOIS - http://addictedtor.free.fr ~~ Etudiant ISUP - CS3 - Industrie et Services ~~http://www.isup.cicrp.jussieu.fr/ ~~ Stagiaire INRIA Futurs - Equipe SELECT ~~ http://www.inria.fr/recherche/equipes/select.fr.html~~ ~ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] CRAN task view for genetics
Hello to everyone! I have built CRAN task view for genetics. For now I have not submit it to CRAN yet and it can be accessible from: http://www.bfro.uni-lj.si/MR/ggorjan/software/R/Genetics.html http://www.bfro.uni-lj.si/MR/ggorjan/software/R/Genetics.ctv I have not submitted it to CRAN, since I would like first some opinion about it. Genetics is really so broad field that I belive one person can not know about every corner. So I would like to get some suggestions, crtitics, pacthes to presentation from my knowledge limited point view. Please note that field of microarray/proteomics, ... area is substantially supported in R and Bioconductor and if someone is prepaired to write a CRAN task view for this I can remove that piece text from my view and provide a link. I have sent this mail also to all mentioned/listed package maintainers. Lep pozdrav / With regards, Gregor Gorjanc -- University of Ljubljana Biotechnical FacultyURI: http://www.bfro.uni-lj.si/MR/ggorjan Zootechnical Department mail: gregor.gorjanc at bfro.uni-lj.si Groblje 3 tel: +386 (0)1 72 17 861 SI-1230 Domzale fax: +386 (0)1 72 17 888 Slovenia, Europe -- One must learn by doing the thing; for though you think you know it, you have no certainty until you try. Sophocles ~ 450 B.C. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Specifying medoids in PAM?
David == David Finlayson [EMAIL PROTECTED] on Tue, 7 Jun 2005 12:11:25 -0700 writes: David I am using the PAM algorithm in the CLUSTER library. David When I allow PAM to seed the medoids using the default __build__ David algorithm things work David well: pam(stats.table, metric=euclidean, stand=TRUE, k=5) David But I have some clusters from a Hierarchical analysis that I would David like to use as seeds for the PAM algorithm. I can't figure what the David mediod argument wants. When I put in the five integer indicies for the David observations in stats.table that I would like to use as seeds (the row David numbers), I segfault R. pam(stats.table, metric=euclidean, stand=TRUE, medoids=c(1,3,20,2,5), k=5) David *** R Crashes *** this is not very helpful. Can you please READ the posting guide and then do as the guide says : - post a *reproducible* example - tell more about what happens when R Crashes David Here is my version info: version David _ David platform i386-pc-mingw32 David arch i386 David os mingw32 David system i386, mingw32 David status David major2 David minor0.1 David year 2004 David month11 David day 15 David language R David Any guidance would be appreciated. David David David -- David David Finlayson David Marine Geology Geophysics David School of Oceanography David Box 357940 David University of Washington David Seattle, WA 98195-7940 David USA David Office: Marine Sciences Building, Room 112 David Phone: (206) 616-9407 David Web: http://students.washington.edu/dfinlays David __ David R-help@stat.math.ethz.ch mailing list David https://stat.ethz.ch/mailman/listinfo/r-help David PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [ESS] Re: [R] Strange characters in 2.1.0?
PaCo == Patrick Connolly [EMAIL PROTECTED]
on Wed, 8 Jun 2005 11:31:44 +1200 writes:
PaCo On Tue, 07-Jun-2005 at 04:10PM +0200, Martin Maechler wrote:
PaCo | Dan == Dan Bolser [EMAIL PROTECTED]
PaCo | ..
PaCo |
PaCo | Dan I have gone back to 2.0.0 :)
PaCo |
PaCo | Don't do that!
PaCo | You've lost tons of nice new features and gained quite an amount
PaCo | of old bugs by downgrading ..
PaCo I get the non-generic quotes to show on the screen, but they won't
PaCo print with enscript. I end up with a lot of wrapped lines and
PaCo nonsense where an unknown character should be.
Why is this diverted from R- to ESS-help?
Printing with enscript is also a topic for printing a
transcript 'output.Rout' resulting e.g. from R CMD BATCH input.R output.Rout
I'm committing a cross-posting felony now, by posting back
to R-help {and please drop ESS-help from cc when further replying}
PaCo What do I need to do to get enscript to know about such characters?
PaCo There is an encoding parameter which defaults to latin1. Should I
PaCo change that to something?
Yes, in principle. latin1 aka ISO-latin-1 aka iso-8859-1
is (for western European languages) the predecessor standard of
the new unicode standard where we use the UTF-8 encoding
{and the above is (too) much simplified; also enter locale
settings and standards}
However, my version of enscript does not seem to support UTF-8 (yet).
Nor does 'a2ps' an alternative to enscript which does pretty
print R source files.
So there are basically two options :
1) Get rid of unicode / utf-8
by setting the locale of your computer / login
to use the old locales, e.g. en_US instead of en_US.utf-8.
This will be more or less fine for Emacs and R --- though in
in our {Redhat Enterprise} setup, the X11-fonts for
non-utf-locales are quite crippled compared to those for
utf-8 ones.
However, as more and more other utilities are based on utf-8
encoded files, you will see funny characters there
if you are using locales like de_* or fr_*, at least,
e.g. for man pages which are only in utf-8 for our Redhat OS setup.
2) Improve the printing tools by
a) filtering *.utf-8 to latin-*
b) printing the resulting latin-*
For filtering, there are programs like 'recode' (was GNU
recode, now Free recode) which are extremely flexible and
'iconv' (less flexible but wider spread) that can translate
utf-8 to and from all kind of encodings / character sets.
In the future, of course everything will work out of the box
when all the utilities in your computer will be aware of utf
encodings and will automatically send correct stuff to the
printer and display it correctly in all kind of viewers/editors... :-)
Given my experiences during the last several months
(where I, e.g., also found that our oldish LaTeX setup
didn't yet accept \usepackage[utf8]{inputencoding ),
If I were in New Zeeland and would not need accents or umlauts,
I'd probably stick with latin1 (and would make sure my X
server got proper non-utf8 fonts) for another year or so.
Martin
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Overlaying Barplots
Hello! I would like to overlay barplot(1:10) with a barplot(seq(1:5, each=2)), indicating that 50% of each bar belongs to category X. How do I do this in R? Best wishes, Sven C. Koehler __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Fitting Theoretical Distributions to Daily Rainfall Data
Dear List Members, I need a bit help about fitting some theoretical distributions (such as geometric, exponential, lognormal or weibull distribution) to the following *dry spell*, *wet spell*, *cycles (Wet-Dry or Dry-Wet)* from my meteorological (daily rainfall) data http://www.angelfire.com/ab5/get5/R.rainfall.txt only for rainy seasen (july - september) of 14 years only: dryspell-c(1,1,4,1,1,4,2,2,3,4,1,2,1,1,1,3,1,3, 2,3,1,2,3,3,2,2,6,2,1,1,3,1,2,1,4,4,1,1,2,1,1,2, 2,4,1,1,1,3,4,1,7,1,3,2,5,1,3,1,1,3,4,8,4,3,1,1, 1,2,3,1,1,2,1,1,2,2,2,3,3,13,13,7,1,1,1,1,7,3,2, 1,3,1,1,2,5,1,2,1,3,1,1,2,3,2,3,2,1,5,1,2,2,2,1, 9,2,2,1,1,4,5,1,1,3,1,3,3,2,1,1,1,17,1,4,5,1,1, 1,1,2,1,2,1,1,7,3,8,2,1,1,2,4,5,1,1,1,2,3,1,1,2, 1,1,3,2,3,1,1,1,3,6,4,1,2,1,2,2,4,2,4,2,1,2,1,3, 1,2,2,1,1,1,1,4,2,1,12,3,1,1,5,1,1,5,1,2,1,1,1,1, 5,3,1,1,3,1,1,6,10,1,1,1,2,1,3,2,2,5,1,1,2,2,1,2, 3,1,3,6,2,1,1,4,6,1,1,1,3,2,2,1,1,1,5,1,1,1,1,3, 1,2,1,7,1,3,1,3,4,1,1,3,4,3,1,4,4,1,3,1,5,3,1,2, 2,1,1,1,2,1,1,6,1,1,1,3,1,3,4,1,1,3,4,1,1,8,1) wetspell-c(1,5,6,4,1,5,3,4,5,2,3,1,5,4,1,4,1,2,3 ,1,5,4,5,2,1,1,1,6,2,19,5,4,6,5,2,7,1,3,1,1,2,1, 3,8,2,3,1,2,5,1,3,8,9,1,1,7,1,2,3,7,9,4,4,1,2,3, 1,1,1,1,1,2,6,7,1,4,1,6,1,5,5,3,2,3,1,1,1,1,6,1, 3,2,1,3,5,6,3,2,6,1,1,3,1,7,3,5,1,2,2,3,1,12,1,8, 3,1,2,1,1,2,1,2,4,2,3,1,1,3,1,4,1,6,5,2,11,6,2,1, 1,9,2,7,1,7,4,1,6,4,8,2,1,1,1,9,3,3,7,2,1,3,3,8,2, 1,7,1,2,2,1,1,1,1,1,5,1,1,3,1,1,1,9,1,7,1,4,3,1,5, 7,1,5,1,5,6,8,5,3,4,1,2,7,9,3,1,4,2,1,1,2,3,1,1,8, 5,2,1,1,1,4,1,1,1,8,4,9,6,3,1,6,5,3,5,2,2,1,5,9,8, 1,6,4,1,2,8,6,1,3,1,2,2,2,3,1,1,5,2,3,11,1,1,1,5, 3,5,1,2,1,9,3,1,1,1,4,10,6,1,1,1,1,1,3,4,1,2,1,5, 2,1,3,2,9,2,1,1,4,2,1,2,9,3,1,1,1,2,6,6,3) cycleWetDry-c(2,9,7,5,5,7,5,7,9,3,5,2,6,5,4,5,4,4, 6,2,7,7,12,5,3,3,7,8,3,20,8,5,8,6,6,11,2,4,3,2,3,3, 5,9,6,4,2,3,8,5,4,15,10,4,3,12,2,5,4,8,12,8,9,5,5, 4,2,2,3,4,2,3,8,8,2,6,3,8,4,8,18,16,9,4,2,2,2,8,9, 3,4,5,2,4,7,11,4,4,7,4,2,4,3,10,5,8,3,3,16,4,3,14, 3,9,12,3,4,2,2,6,6,3,5,5,4,4,4,5,2,5,2,17,6,6,16, 7,3,2,2,11,3,9,2,8,11,4,14,6,16,3,3,5,6,10,4,4,9, 5,2,4,5,9,3,4,9,4,3,3,2,4,7,5,2,10,2,3,5,5,3,5,11, 2,9,2,7,4,3,7,8,2,6,2,9,8,9,9,14,6,5,2,7,8,10,8,2, 6,3,2,2,3,8,4,2,9,8,3,2,7,3,5,2,2,3,9,7,11,8,8,2, 7,7,5,6,4,5,2,8,12,10,2,7,8,7,3,9,7,4,5,3,3,3,3,8, 2,2,6,3,6,18,2,3,2,12,4,8,2,5,5,10,4,4,5,4,5,14, 10,2,4,2,6,4,4,6,3,4,2,6,3,3,4,3,15,3,2,2,7,3,4, 6,10,4,4,5,2,3,14,7,4) cycleDryWet-c(2,6,10,5,2,9,5,6,8,6,4,3,6,5,2,7,2, 5,5,4,6,6,8,4,4,3,3,12,4,20,6,7,7,7,3,11,5,4,2,3,3, 2,5,10,3,7,2,3,6,4,7,9,16,2,4,9,6,3,6,8,10,7,8,9,6, 6,2,2,2,3,4,3,7,9,2,5,3,8,3,8,8,16,15,11,2,2,2,2,13, 4,5,3,4,4,6,8,8,3,8,2,4,4,2,9,6,7,4,4,3,8,2,14,3,10, 4,10,4,3,2,3,5,7,5,3,6,2,4,6,3,5,2,7,22,3,15,11,3,2, 2,10,4,8,3,8,5,8,9,12,10,3,2,3,5,14,4,4,8,4,4,4,4, 10,3,2,10,3,5,3,2,2,4,7,5,6,3,2,5,3,5,3,13,3,8,3, 5,6,2,7,9,2,6,2,6,10,10,6,3,15,7,2,3,12,10,4,6,5, 4,2,2,3,4,6,4,9,6,5,2,2,7,14,2,2,2,10,5,12,8,5,6, 7,6,5,7,3,4,4,6,12,14,3,7,5,5,8,9,7,2,6,3,4,3,3,4, 6,2,6,3,4,14,15,2,3,6,10,6,4,3,4,13,4,2,4,5,7,11, 10,5,2,4,2,6,6,5,3,4,3,6,3,2,5,3,10,8,2,2,5,5,2, 5,13,4,2,4,5,3,7,14,4) Using table() to each dryspell, wetspell, cycleWetDry, cycleDryWet we find the empirical distribution functions all of which seem to be positively skewed with long tail. Therefore, i'd like to fit geometric, exponential, lognormal or weibull distribution for each dryspell, wetspell, cycleWetDry, cycleDryWet. Better fit may be defined by higher p-values of goodness-of-fit tests. Is there any way i can do fit data to those theoretical distributions in R? Is there any existing program/function/package to solve such problem? Any suggestion, direction, references, help, replies will be highly appreciated. Thank you for your time. -- Mohammad Ehsanul Karim Web: http://snipurl.com/ehsan ISRT, University of Dhaka, BD -- __ Get on-the-go sports scores, stock quotes, news and more. Check it out! __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Fitting Theoretical Distributions to Daily Rainfall Data
Have a look at the fit.dist function in Jim Lindsey's gnlm package at
http://popgen0146uns50.unimaas.nl/~jlindsey/rcode.html
fit.dist {gnlm} R Documentation
Fit Probability Distributions to Frequency Data
Description
fit.dist fits the distributions in Chapter 4 of Lindsey (1995, 2003 2nd edn):
binomial, beta-binomial, Poisson, negative binomial, geometric, zeta, normal,
log normal, inverse Gauss, logistic, Laplace, Cauchy, Student t, exponential,
Pareto, gamma, and Weibull to frequency (histogram) data, possibly plotting
the frequency polygon of fitted values with the histogram.
fitdistr from the MASS package works quite well, too.
Dear List Members,
I need a bit help about fitting some theoretical
distributions (such as geometric, exponential,
lognormal or weibull distribution)
Ken Knoblauch
Inserm U371, Cerveau et Vision
Department of Cognitive Neurosciences
18 avenue du Doyen Lepine
69500 Bron
France
tel: +33 (0)4 72 91 34 77
fax: +33 (0)4 72 91 34 61
portable: 06 84 10 64 10
http://www.lyon.inserm.fr/371/
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] bug in predict.lme?
Dear All, I've come across a problem in predict.lme. Assigning a model formula to a variable and then using this variable in lme (instead of typing the formula into the formula part of lme) works as expect. However, when performing a predict on the fitted model I gan an error messag - predict.lme (but not predictlm) seems to expect a 'properly' typed in formula and a cannot extract the formula from the variable. THe code below demonstrates this. Is this a known or expected behavour of predict.lme or is this a bug? kind regards, Arne (R-2.1.0) library(nlme) ... mod - distance ~ age + Sex # example from ?lme mod distance ~ age + Sex fm2 - lme(mod, data = Orthodont, random = ~ 1) anova(fm2) numDF denDF F-value p-value (Intercept) 180 4123.156 .0001 age 180 114.838 .0001 Sex 1259.292 0.0054 fm2 Linear mixed-effects model fit by REML Data: Orthodont Log-restricted-likelihood: -218.7563 Fixed: mod ... predict(fm2, Orthodont) Error in mCall[[fixed]][-2] : object is not subsettable fm2 - update(fm2, . ~ .) # this replaces mod by the contents of variable mod fm2 Linear mixed-effects model fit by REML Data: Orthodont Log-restricted-likelihood: -218.7563 Fixed: distance ~ age + Sex ... predict(fm2, Orthodont) M01 M01 M01 M01 ... 25.39237 26.71274 28.03311 29.35348 21.61052 ... fm2 - lm(mod, data = Orthodont) predict(fm2, Orthodont) 1234 ... 22.98819 24.30856 25.62894 26.94931 ... __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] How to change the value of a class slot
Ross Boylan wrote:
I defined an S4 class with a slot i. Then I wrote a regular function
that attempted to increment i.
[... details deleted ...]
What do I need to do to update slot values?
Here are some possibly relevant code fragments
setClass(CompletePathMaker,
representation(i=integer,
timeOffset=numeric, # to avoid 0's
truePaths=TruePaths)
)
nextPath - function(pm){ #pm is a CompletePathMaker
[EMAIL PROTECTED] - [EMAIL PROTECTED](1)
[etc]
I'm trying to make the class behave like an iterator, with i keeping
track of its location. I'm sure there are more R'ish ways to go, but
I'm also pretty sure I'm going to want to be able to update slots.
Hello Ross,
I see that your question was related to S4, but I just noticed a
solution based on the R.oo package so I thought I would add a solution
based on the proto package too. We had similar problems several times
ago and (to my surprise) found R to be an extremely flexible language
even for these things. Our favorite solution is available as
package(proto). It requires R 2.1, because of several subtle
improvements regarding environments, which made our implementation more
streamlined.
Does the following example do what you intended?
##=
library(proto)
## 1) define an object
CompletePathMaker - proto(
index = 0,
bumpIndex = function(., dindex = 1)
.$index - .$index + as.integer(dindex)
)
## 2) create a child object of CompletePathMaker
cpm - CompletePathMaker$proto()
## 3) set the index component to 3
cpm$index - 3
## 4) iterate the index
cpm$bumpIndex(2)
## print the result
cpm$index
##=
This approach is very compact and needs only one new function: proto.
Also note how simple it is conceptually. We did not even create any
classes. We just created a parent object CompletePathMaker and a child
to it, cpm, and got everything else via delegation (i.e. inheritance).
Hope it helps
Thomas P.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Robustness of Segmented Regression Contributed by Muggeo
Hello, R users, I applied segmented regression method contributed by Muggeo and got different slope estimates depending on the initial break points. The results are listed below and I'd like to know what is a reasonable approach handling this kinds of problem. I think applying various initial break points is certainly not a efficient approach. Is there any other methods to deal with segmented regression? From a graph, v shapes are more clear at 1.2 and 1.5 break points than 1.5 and 1.7. Appreciate your help. Result1: Initial break points are 1.2 and 1.5. The estimated break points and slopes: Estimated Break-Point(s): Est. St.Err Mean.Vel 1.285 0.05258 1.6520.01247 Est. St.Err. t valueCI(95%).l CI(95%).u slope1 0.4248705 0.3027957 1.403159-0.16859821.018339 slope2 2.3281445 0.3079903 7.559149 1.72449462.931794 slope3 9.5425516 0.7554035 12.632390 8.0619879 11.023115 Adjusted R-squared: 0.9924. Result2: Initial break points are 1.5 and 1.7. The estimated break points and slopes: Estimated Break-Point(s): Est. St.Err Mean.Vel 1.412 0.02195 1.699 0.01001 Est. St.Err.t valueCI(95%).l CI(95%).u slope1 0.7300483 0.13815875.284129 0.4592623 1.000834 slope2 3.4479466 0.244253014.116289 2.9692194 3.926674 slope3 12.500 1.7783840 7.028853 9.0144314 15.985569 Adjusted R-squared: 0.995. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Overlaying Barplots
On Wed, 2005-06-08 at 12:07 +0200, Sven C. Koehler wrote: Hello! I would like to overlay barplot(1:10) with a barplot(seq(1:5, each=2)), indicating that 50% of each bar belongs to category X. How do I do this in R? If you pass a matrix to barplot, it will stack values from the same column on top of each other. So something like barplot(rbind(1:10, 1:10)) will do what you want. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Robustness of Segmented Regression Contributed by Muggeo
You might try rqss() in the quantreg package. It gives piecewise linear fits for a nonparametric form of median regression using total variation of the derivative of the fitted function as a penalty term. A tuning parameter (lambda) controls the number of distinct segments. More details are available in the vignette for the quantreg package. url:www.econ.uiuc.edu/~rogerRoger Koenker email[EMAIL PROTECTED]Department of Economics vox: 217-333-4558University of Illinois fax: 217-244-6678Champaign, IL 61820 On Jun 8, 2005, at 7:25 AM, Park, Kyong H Mr. RDECOM wrote: Hello, R users, I applied segmented regression method contributed by Muggeo and got different slope estimates depending on the initial break points. The results are listed below and I'd like to know what is a reasonable approach handling this kinds of problem. I think applying various initial break points is certainly not a efficient approach. Is there any other methods to deal with segmented regression? From a graph, v shapes are more clear at 1.2 and 1.5 break points than 1.5 and 1.7. Appreciate your help. Result1: Initial break points are 1.2 and 1.5. The estimated break points and slopes: Estimated Break-Point(s): Est. St.Err Mean.Vel 1.285 0.05258 1.6520.01247 Est. St.Err. t valueCI (95%).l CI(95%).u slope1 0.4248705 0.3027957 1.403159-0.1685982 1.018339 slope2 2.3281445 0.3079903 7.559149 1.7244946 2.931794 slope3 9.5425516 0.7554035 12.632390 8.0619879 11.023115 Adjusted R-squared: 0.9924. Result2: Initial break points are 1.5 and 1.7. The estimated break points and slopes: Estimated Break-Point(s): Est. St.Err Mean.Vel 1.412 0.02195 1.699 0.01001 Est. St.Err.t valueCI (95%).l CI(95%).u slope1 0.7300483 0.13815875.284129 0.4592623 1.000834 slope2 3.4479466 0.244253014.116289 2.9692194 3.926674 slope3 12.500 1.7783840 7.028853 9.0144314 15.985569 Adjusted R-squared: 0.995. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting- guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Robustness of Segmented Regression Contributed by Muggeo
On Wed, 8 Jun 2005 08:25:16 -0400 Park, Kyong H Mr. RDECOM wrote: Hello, R users, I applied segmented regression method contributed by Muggeo and got different slope estimates depending on the initial break points. The results are listed below and I'd like to know what is a reasonable approach handling this kinds of problem. I think applying various initial break points is certainly not a efficient approach. Is there any other methods to deal with segmented regression? From a graph, v shapes are more clear at 1.2 and 1.5 break points than 1.5 and 1.7. Appreciate your help. When you keep the number of break points fixed, then there is a unique solution to the problem of fitting a segmented regression: the solution which maximizes the likelihood (or for linear models equivalently minimizes the RSS). Vito's segmented package gives an iterative method which can be shown to converge to this unique solution. If empirically you find different solutions with different starting values, you can always compare them using the RSS or log-likelihood and choose the one which fits better (because the other one can't be the optimal solution). The function breakpoints() in package strucchange computes (as opposed to approximates) the unique solution for a fully segmented model instead of a broken line trend. Another nonparametric solution using quantreg was already pointed out by Roger. hth, Z Result1: Initial break points are 1.2 and 1.5. The estimated break points and slopes: Estimated Break-Point(s): Est. St.Err Mean.Vel 1.285 0.05258 1.6520.01247 Est. St.Err. t value CI(95%).l CI(95%).u slope1 0.4248705 0.3027957 1.403159-0.1685982 1.018339 slope2 2.3281445 0.3079903 7.559149 1.7244946 2.931794 slope3 9.5425516 0.7554035 12.632390 8.0619879 11.023115 Adjusted R-squared: 0.9924. Result2: Initial break points are 1.5 and 1.7. The estimated break points and slopes: Estimated Break-Point(s): Est. St.Err Mean.Vel 1.412 0.02195 1.699 0.01001 Est. St.Err.t value CI(95%).l CI(95%).u slope1 0.7300483 0.13815875.284129 0.4592623 1.000834 slope2 3.4479466 0.244253014.116289 2.9692194 3.926674 slope3 12.500 1.7783840 7.028853 9.0144314 15.985569 Adjusted R-squared: 0.995. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Random seed problem in MCMC coupling of chains
Hello!
I am performing coupling of chains in MCMC and I need the same value
of seed for two chains. I will show demo of what I want:
R code, which might show my example is:
niter - 3
nchain - 2
tmpSeed - 123
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(tmpSeed)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
tmpSeed - .Random.seed
}
}
I get this:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.43588
iter: 2 chain: 2 runif: 0.43588
iter: 3 chain: 1 runif: 0.43588
iter: 3 chain: 2 runif: 0.43588
but I would like to get:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.67676
iter: 2 chain: 2 runif: 0.67676
iter: 3 chain: 1 runif: 0.12368
iter: 3 chain: 2 runif: 0.12368
Note that seed value is of course changing, but it is parallel
between chains.
I am able to do only this, since setting seed at the beginning
of chain i.e iteration is not a problem, but I want an upper
scheme, since I compare chains and stop one if some condition is
satisfied.
tmpSeed - 123
for (i in 1:nchain) { # chains
set.seed(tmpSeed)
for (j in 1:niter) { # iterations
a - runif(1)
cat(iter:, j, chain:, i, runif:, a, \n)
}
}
iter: 1 chain: 1 runif: 0.28758
iter: 2 chain: 1 runif: 0.7883
iter: 3 chain: 1 runif: 0.40898
iter: 1 chain: 2 runif: 0.28758
iter: 2 chain: 2 runif: 0.7883
iter: 3 chain: 2 runif: 0.40898
iter: 1 chain: 3 runif: 0.28758
iter: 2 chain: 3 runif: 0.7883
iter: 3 chain: 3 runif: 0.40898
I was looking in 'rlecuyer', 'rsprng' and 'setRNG', but did not find
anything usable for me. From reading on http://sprng.cs.fsu.edu/
'rsprng' provides just opposite of what I want, 'rlecuyer' is a bit
to technical for me, but I think it also doesn't give identical
seed for parallels. 'setRNG', especially it's function 'getRNG()'
looks nice but its arguments should have seed stored. How can one
do that?
Thanks in advance!
Lep pozdrav / With regards,
Gregor Gorjanc
--
University of Ljubljana
Biotechnical FacultyURI: http://www.bfro.uni-lj.si/MR/ggorjan
Zootechnical Department mail: gregor.gorjanc at bfro.uni-lj.si
Groblje 3 tel: +386 (0)1 72 17 861
SI-1230 Domzale fax: +386 (0)1 72 17 888
Slovenia, Europe
--
One must learn by doing the thing; for though you think you know it,
you have no certainty until you try. Sophocles ~ 450 B.C.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Converting code from MATLAB to R
Hi, I'm having trouble converting code from MATLAB to R; I want to find the equivalence to MATLAB's function 'fsolve'. I've tried 'nlm', on the squared argument, in R but i did not get the same results. Thankful if helped. Best regards, Martin Englund -- This e-mail and any attachment may be confidential and may also be privileged. If you are not the intended recipient, please notify us immediately and then delete this e-mail and any attachment without retaining copies or disclosing the contents thereof to any other person. Thank you. -- [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] RODBC question: problem importing series with blank cells
Hello,
I have an excel file that I load through RODBC. Some of my columns are
blank. They are equity time series and the stocks did not exist at the
earlier dates. My problem is that the whole column becomes NA even though
there are numbers at later dates.
Here's my excel file
http://www.tradebit.com/download.php/35699
http://www.tradebit.com/download.php/35699
And here's the code I use:
library(RODBC)
chan - odbcConnectExcel(C:/book54.xls) #load data
ts- sqlFetch(chan, Sheet1)
close(chan)
ts-ts[-1,]
str(ts)
head(ts)
tail(ts)
Regards,
Pierre Lapointe
Assistant Market Strategist
National Bank Financial
***
AVIS DE NON-RESPONSABILITE:\ Ce document transmis par courri...{{dropped}}
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Random seed problem in MCMC coupling of chains
do you want something like this:
niter - 3
nchain - 2
rs - sample(500, niter, TRUE)
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(rs[i])
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
}
}
I hope it helps.
Best,
Dimitris
Dimitris Rizopoulos
Ph.D. Student
Biostatistical Centre
School of Public Health
Catholic University of Leuven
Address: Kapucijnenvoer 35, Leuven, Belgium
Tel: +32/16/336899
Fax: +32/16/337015
Web: http://www.med.kuleuven.ac.be/biostat/
http://www.student.kuleuven.ac.be/~m0390867/dimitris.htm
- Original Message -
From: Gorjanc Gregor [EMAIL PROTECTED]
To: r-help@stat.math.ethz.ch
Sent: Wednesday, June 08, 2005 3:27 PM
Subject: [R] Random seed problem in MCMC coupling of chains
Hello!
I am performing coupling of chains in MCMC and I need the same value
of seed for two chains. I will show demo of what I want:
R code, which might show my example is:
niter - 3
nchain - 2
tmpSeed - 123
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(tmpSeed)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
tmpSeed - .Random.seed
}
}
I get this:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.43588
iter: 2 chain: 2 runif: 0.43588
iter: 3 chain: 1 runif: 0.43588
iter: 3 chain: 2 runif: 0.43588
but I would like to get:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.67676
iter: 2 chain: 2 runif: 0.67676
iter: 3 chain: 1 runif: 0.12368
iter: 3 chain: 2 runif: 0.12368
Note that seed value is of course changing, but it is parallel
between chains.
I am able to do only this, since setting seed at the beginning
of chain i.e iteration is not a problem, but I want an upper
scheme, since I compare chains and stop one if some condition is
satisfied.
tmpSeed - 123
for (i in 1:nchain) { # chains
set.seed(tmpSeed)
for (j in 1:niter) { # iterations
a - runif(1)
cat(iter:, j, chain:, i, runif:, a, \n)
}
}
iter: 1 chain: 1 runif: 0.28758
iter: 2 chain: 1 runif: 0.7883
iter: 3 chain: 1 runif: 0.40898
iter: 1 chain: 2 runif: 0.28758
iter: 2 chain: 2 runif: 0.7883
iter: 3 chain: 2 runif: 0.40898
iter: 1 chain: 3 runif: 0.28758
iter: 2 chain: 3 runif: 0.7883
iter: 3 chain: 3 runif: 0.40898
I was looking in 'rlecuyer', 'rsprng' and 'setRNG', but did not find
anything usable for me. From reading on http://sprng.cs.fsu.edu/
'rsprng' provides just opposite of what I want, 'rlecuyer' is a bit
to technical for me, but I think it also doesn't give identical
seed for parallels. 'setRNG', especially it's function 'getRNG()'
looks nice but its arguments should have seed stored. How can one
do that?
Thanks in advance!
Lep pozdrav / With regards,
Gregor Gorjanc
--
University of Ljubljana
Biotechnical FacultyURI:
http://www.bfro.uni-lj.si/MR/ggorjan
Zootechnical Department mail: gregor.gorjanc at bfro.uni-lj.si
Groblje 3 tel: +386 (0)1 72 17 861
SI-1230 Domzale fax: +386 (0)1 72 17 888
Slovenia, Europe
--
One must learn by doing the thing; for though you think you know
it,
you have no certainty until you try. Sophocles ~ 450 B.C.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Random seed problem in MCMC coupling of chains
On 6/8/2005 9:27 AM, Gorjanc Gregor wrote:
Hello!
I am performing coupling of chains in MCMC and I need the same value
of seed for two chains. I will show demo of what I want:
R code, which might show my example is:
niter - 3
nchain - 2
tmpSeed - 123
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(tmpSeed)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
tmpSeed - .Random.seed
}
}
I get this:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.43588
iter: 2 chain: 2 runif: 0.43588
iter: 3 chain: 1 runif: 0.43588
iter: 3 chain: 2 runif: 0.43588
but I would like to get:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.67676
iter: 2 chain: 2 runif: 0.67676
iter: 3 chain: 1 runif: 0.12368
iter: 3 chain: 2 runif: 0.12368
Note that seed value is of course changing, but it is parallel
between chains.
set.seed takes an integer, but .Random.seed is a complicated vector.
You need to play with .Random.seed directly, and move your setting of
tmpSeed out of the inner loop, i.e.
niter - 3
nchain - 2
set.seed(123)
tmpSeed - .Random.seed
for (i in 1:niter) { # iterations
+ for (j in 1:nchain) { # chains
+ .Random.seed - tmpSeed
+ a - runif(1)
+ cat(iter:, i, chain:, j, runif:, a, \n)
+ }
+ tmpSeed - .Random.seed
+ }
iter: 1 chain: 1 runif: 0.2875775
iter: 1 chain: 2 runif: 0.2875775
iter: 2 chain: 1 runif: 0.7883051
iter: 2 chain: 2 runif: 0.7883051
iter: 3 chain: 1 runif: 0.4089769
iter: 3 chain: 2 runif: 0.4089769
However, heed the warnings in ?set.seed: with some generators
.Random.seed *does not* contain the full state of the RNG. As far as I
know there is no way to obtain the full state.
Duncan Murdoch
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Converting code from MATLAB to R
Hi, I'm having trouble converting code from MATLAB to R; I want to find the equivalence to MATLAB's function 'fsolve'. I've tried 'nlm', on the squared argument, in R but i did not get the same results. Thankful if helped. Best regards, Martin Englund -- This e-mail and any attachment may be confidential and may also be privileged. If you are not the intended recipient, please notify us immediately and then delete this e-mail and any attachment without retaining copies or disclosing the contents thereof to any other person. Thank you. -- [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to run a script at the beginning of an interacive session ?
On Tue, 2005-06-07 at 21:33 -0700, Mike R wrote:
On 6/7/05, Marc Schwartz [EMAIL PROTECTED] wrote:
On 6/7/05, Gabor Grothendieck [EMAIL PROTECTED] wrote:
On 6/7/05, Liaw, Andy [EMAIL PROTECTED] wrote:
snips
Thanks Andy, Gabor and Marc.
contents of .Rprofile
.First - function() {source(startup.r)}
contents of wrapR.v1
#!/bin/bash
ln -s $1 startup.r
R
rm -f startup.r
contents of wrapR.v2
#!/bin/bash
echo .First - function() {source(\$1\)} .Rprofile
R
contents of project_A.r
library(graphics)
plot(1:10)
and then start R with either of the following command lines:
## ./wrapR.v1 project_A.r
## ./wrapR.v2 project_A.r
What do you think? is this the most simple way?
Thanks in advance,
Mike
Mike,
There are several options here, and I would try to stay with a simple
approach, minimizing having to modify files, such as you are doing in
wrapR.v2, which BTW will overwrite all other contents of ~/.Rprofile.
You could get very complicated if you wanted and use something like
'sed' to engage in non-interactive editing of the ~./Rprofile file, but
I think that would increase the complexity beyond what is needed here.
One possible approach, subject to comment by others, would be to have
your normal default ~/.Rprofile file containing your base specifications
for a default startup.
Then, create one or more ~/.RprofileX files, where X is a number of 1:n,
for each separate file.
In each of these files, fully construct a .First function (as opposed to
sourcing() a separate file), such as:
# In file ~/.Rprofile1
...
.First - function()
{
library(graphics)
plot(1:10)
}
# In file ~/.Rprofile2
...
.First - function()
{
library(graphics)
barplot(1:10)
}
In each of the above, replace the '...' with the core contents of
~/.Rprofile so that besides .First, all else is the same if you want it
that way.
With the above in place, now run R from a terminal or within a shell
script using:
R_PROFILE=~/.Rprofile1 R --no-init-file
OR
R_PROFILE=~/.Rprofile2 R --no-init-file
In this way, by setting the R_PROFILE environment variable, you can
define which init file to use on startup. The use of --no-init-file
ignores your default ~/.Rprofile file to avoid any possible conflicts
there.
Note, importantly, this approach would obviate the use of any site-
wide .Rprofile file, which by default would be
$R_HOME/etc/Rprofile.site.
Now all you have to do is to create a new ~/.RprofileX each time you
want to utilize a different .First by modifying the command line call to
R.
If you want to run R with a default startup (ie. no .First function),
just use R as per normal.
HTH,
Marc Schwartz
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] RODBC question: problem importing series with blank cells
On Wed, 8 Jun 2005, Lapointe, Pierre wrote:
Hello,
I have an excel file that I load through RODBC. Some of my columns are
blank. They are equity time series and the stocks did not exist at the
earlier dates. My problem is that the whole column becomes NA even though
there are numbers at later dates.
My Excel ODBC driver is sending NULLs for those rows, so I think this is
an ODBC driver problem.
Here's my excel file
http://www.tradebit.com/download.php/35699
http://www.tradebit.com/download.php/35699
And here's the code I use:
library(RODBC)
chan - odbcConnectExcel(C:/book54.xls) #load data
ts- sqlFetch(chan, Sheet1)
close(chan)
ts-ts[-1,]
str(ts)
head(ts)
tail(ts)
Regards,
Pierre Lapointe
Assistant Market Strategist
National Bank Financial
***
AVIS DE NON-RESPONSABILITE:\ Ce document transmis par courri...{{dropped}}
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
--
Brian D. Ripley, [EMAIL PROTECTED]
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@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Converting code from MATLAB to R
Hi Martin optim() should do the trick [and this isn't in R-and-octave2.txt; I'll update when I get a minute] HTH rksh On Jun 8, 2005, at 02:27 pm, Martin Englund wrote: Hi, I'm having trouble converting code from MATLAB to R; I want to find the equivalence to MATLAB's function 'fsolve'. I've tried 'nlm', on the squared argument, in R but i did not get the same results. Thankful if helped. Best regards, Martin Englund --- --- This e-mail and any attachment may be confidential and may also be privileged. If you are not the intended recipient, please notify us immediately and then delete this e-mail and any attachment without retaining copies or disclosing the contents thereof to any other person. Thank you. --- --- [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Robin Hankin Uncertainty Analyst National Oceanography Centre, Southampton European Way, Southampton SO14 3ZH, UK tel 023-8059-7743 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Bounding or constraining parameters in non-linear regressions
On 6/7/05, [EMAIL PROTECTED] [EMAIL PROTECTED] wrote: Dear R-Users, Being an engineer and not a statistician, my desired course of action may either be impossible or very simple. I am attempting to fit a non-linear model to some measured data. One term in the model contains a square-root, but in the course of regression, this term turns negative and an error occurs. I started using Micrsoft's Excel Solver, and then I turned to NIST's Datplot statistical package. I can constrain in Solver, but it violates those constraints. :) Dataplot does not have the capability to constrain parameters. Does R have the capability to constrain or bound parameters in non-linear regressions? Sort of. If you look at the stats package in r-devel you will see that a function called nlminb has been added. This function calls optimization software from the Port package (http://www.netlib.com/port/). The Fortran code for constrained nonlinear least squares problems is included in the package but the interface code for R has not yet been written. The energetic could create such interface code by emulating that for nlminb - it's not that long. Alternatively you could use either optim or nlminb on the function which is the residual sum of squares from your model. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] variable png-size in multiple figure environment
If I put only one plot in a png, I can drag the generated png on the edges and it changes its size by this. But when I put several plots on the same sheet using mfrow, the size can no longer be changed when viewing the file, and the resolution is bad. What do I need to do to keep the variability of size in a multiple figure environment? __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] get level combinations from by list
Dear useRs,
Given this code I end up with a list of class by:
a - sample(1:5,200,replace=TRUE)
b - sample(c(v1,v2,v3),200,replace=TRUE)
c - sample(c(11,22,33),200,replace=TRUE)
data - runif(200)
grouped - by(data,list(a,b,c),function(x) {c(min=min(x),max=max(x),
median=round(median(x),digits=2),mean=round(mean(x),digits=2))})
dfr - do.call(rbind,grouped)#the levels are missing
#--
grouped
typeof(grouped)
class(grouped)
dimnames(grouped)
How do I get at the levels of the 'group by' variables for each
subset? For example, from this part of the by list I want 4, v2 and
33:
: 4
: v2
: 33
min maxmedian mean
0.3897450 0.9215315 0.730 0.670
---
Thank you,
b.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] Random seed problem in MCMC coupling of chains
Thanks to Duncan, Dimitris as well as James for answers. I'll provide
here also example from James, which seems to be the easiest of them
all and was not posted to the list:
niter - 3
nchain - 2
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(i)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
}
}
Note that seed is set with iteration counter. This is really neat and
simple. I am just wondering if this is OK from RNG point of view. Can
someone comment on that?
Lep pozdrav / With regards,
Gregor Gorjanc
--
University of Ljubljana
Biotechnical FacultyURI: http://www.bfro.uni-lj.si/MR/ggorjan
Zootechnical Department mail: gregor.gorjanc at bfro.uni-lj.si
Groblje 3 tel: +386 (0)1 72 17 861
SI-1230 Domzale fax: +386 (0)1 72 17 888
Slovenia, Europe
--
One must learn by doing the thing; for though you think you know it,
you have no certainty until you try. Sophocles ~ 450 B.C.
--
-Original Message-
From: Duncan Murdoch [mailto:[EMAIL PROTECTED]
Sent: sre 2005-06-08 15:53
To: Gorjanc Gregor
Cc: r-help@stat.math.ethz.ch
Subject: Re: [R] Random seed problem in MCMC coupling of chains
On 6/8/2005 9:27 AM, Gorjanc Gregor wrote:
Hello!
I am performing coupling of chains in MCMC and I need the same value
of seed for two chains. I will show demo of what I want:
R code, which might show my example is:
niter - 3
nchain - 2
tmpSeed - 123
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(tmpSeed)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
tmpSeed - .Random.seed
}
}
I get this:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.43588
iter: 2 chain: 2 runif: 0.43588
iter: 3 chain: 1 runif: 0.43588
iter: 3 chain: 2 runif: 0.43588
but I would like to get:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.67676
iter: 2 chain: 2 runif: 0.67676
iter: 3 chain: 1 runif: 0.12368
iter: 3 chain: 2 runif: 0.12368
Note that seed value is of course changing, but it is parallel
between chains.
set.seed takes an integer, but .Random.seed is a complicated vector.
You need to play with .Random.seed directly, and move your setting of
tmpSeed out of the inner loop, i.e.
niter - 3
nchain - 2
set.seed(123)
tmpSeed - .Random.seed
for (i in 1:niter) { # iterations
+ for (j in 1:nchain) { # chains
+ .Random.seed - tmpSeed
+ a - runif(1)
+ cat(iter:, i, chain:, j, runif:, a, \n)
+ }
+ tmpSeed - .Random.seed
+ }
iter: 1 chain: 1 runif: 0.2875775
iter: 1 chain: 2 runif: 0.2875775
iter: 2 chain: 1 runif: 0.7883051
iter: 2 chain: 2 runif: 0.7883051
iter: 3 chain: 1 runif: 0.4089769
iter: 3 chain: 2 runif: 0.4089769
However, heed the warnings in ?set.seed: with some generators
.Random.seed *does not* contain the full state of the RNG. As far as I
know there is no way to obtain the full state.
Duncan Murdoch
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Random seed problem in MCMC coupling of chains
The tools in setRNG are intended for this kind of problem and I do use
them regularly in much more complicated situations. They help save all
the information, in addition to the seed, that you need for reproducible
simulations. Try
niter - 3
nchain - 2
for (i in 1:niter) { # iterations
tmpSeed - setRNG()
for (j in 1:nchain) { # chains
setRNG(tmpSeed)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
}
}
iter: 1 chain: 1 runif: 0.8160078
iter: 1 chain: 2 runif: 0.8160078
iter: 2 chain: 1 runif: 0.4909793
iter: 2 chain: 2 runif: 0.4909793
iter: 3 chain: 1 runif: 0.4425924
iter: 3 chain: 2 runif: 0.4425924
HTH,
Paul Gilbert
Gorjanc Gregor wrote:
Hello!
I am performing coupling of chains in MCMC and I need the same value
of seed for two chains. I will show demo of what I want:
R code, which might show my example is:
niter - 3
nchain - 2
tmpSeed - 123
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(tmpSeed)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
tmpSeed - .Random.seed
}
}
I get this:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.43588
iter: 2 chain: 2 runif: 0.43588
iter: 3 chain: 1 runif: 0.43588
iter: 3 chain: 2 runif: 0.43588
but I would like to get:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.67676
iter: 2 chain: 2 runif: 0.67676
iter: 3 chain: 1 runif: 0.12368
iter: 3 chain: 2 runif: 0.12368
Note that seed value is of course changing, but it is parallel
between chains.
I am able to do only this, since setting seed at the beginning
of chain i.e iteration is not a problem, but I want an upper
scheme, since I compare chains and stop one if some condition is
satisfied.
tmpSeed - 123
for (i in 1:nchain) { # chains
set.seed(tmpSeed)
for (j in 1:niter) { # iterations
a - runif(1)
cat(iter:, j, chain:, i, runif:, a, \n)
}
}
iter: 1 chain: 1 runif: 0.28758
iter: 2 chain: 1 runif: 0.7883
iter: 3 chain: 1 runif: 0.40898
iter: 1 chain: 2 runif: 0.28758
iter: 2 chain: 2 runif: 0.7883
iter: 3 chain: 2 runif: 0.40898
iter: 1 chain: 3 runif: 0.28758
iter: 2 chain: 3 runif: 0.7883
iter: 3 chain: 3 runif: 0.40898
I was looking in 'rlecuyer', 'rsprng' and 'setRNG', but did not find
anything usable for me. From reading on http://sprng.cs.fsu.edu/
'rsprng' provides just opposite of what I want, 'rlecuyer' is a bit
to technical for me, but I think it also doesn't give identical
seed for parallels. 'setRNG', especially it's function 'getRNG()'
looks nice but its arguments should have seed stored. How can one
do that?
Thanks in advance!
Lep pozdrav / With regards,
Gregor Gorjanc
--
University of Ljubljana
Biotechnical FacultyURI: http://www.bfro.uni-lj.si/MR/ggorjan
Zootechnical Department mail: gregor.gorjanc at bfro.uni-lj.si
Groblje 3 tel: +386 (0)1 72 17 861
SI-1230 Domzale fax: +386 (0)1 72 17 888
Slovenia, Europe
--
One must learn by doing the thing; for though you think you know it,
you have no certainty until you try. Sophocles ~ 450 B.C.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] get level combinations from by list
Here is a different approach to achieving what I think you want using
summarize() in the Hmisc package:
library(Hmisc)
mydata - data.frame(a = sample(1:5,200,replace=TRUE),
b = sample(c(v1,v2,v3),200,replace=TRUE),
c = sample(c(11,22,33),200,replace=TRUE),
y = runif(200))
attach(mydata)
summarize(y, llist(a, b, c),
function(x){c(min=min(x),
max=max(x),
median=round(median(x),digits=2),
mean=round(mean(x),digits=2))},
stat.name=min)
detach(mydata)
bogdan romocea wrote:
Dear useRs,
Given this code I end up with a list of class by:
a - sample(1:5,200,replace=TRUE)
b - sample(c(v1,v2,v3),200,replace=TRUE)
c - sample(c(11,22,33),200,replace=TRUE)
data - runif(200)
grouped - by(data,list(a,b,c),function(x) {c(min=min(x),max=max(x),
median=round(median(x),digits=2),mean=round(mean(x),digits=2))})
dfr - do.call(rbind,grouped)#the levels are missing
#--
grouped
typeof(grouped)
class(grouped)
dimnames(grouped)
How do I get at the levels of the 'group by' variables for each
subset? For example, from this part of the by list I want 4, v2 and
33:
: 4
: v2
: 33
min maxmedian mean
0.3897450 0.9215315 0.730 0.670
---
Thank you,
b.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 452-1424 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Random seed problem in MCMC coupling of chains
Beware that your easy trick will give you the same result every time you
run it. You need a better scheme if you actually intend to get a new
experiment each time you run it.
Paul
Gorjanc Gregor wrote:
Thanks to Duncan, Dimitris as well as James for answers. I'll provide
here also example from James, which seems to be the easiest of them
all and was not posted to the list:
niter - 3
nchain - 2
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(i)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
}
}
Note that seed is set with iteration counter. This is really neat and
simple. I am just wondering if this is OK from RNG point of view. Can
someone comment on that?
Lep pozdrav / With regards,
Gregor Gorjanc
--
University of Ljubljana
Biotechnical FacultyURI: http://www.bfro.uni-lj.si/MR/ggorjan
Zootechnical Department mail: gregor.gorjanc at bfro.uni-lj.si
Groblje 3 tel: +386 (0)1 72 17 861
SI-1230 Domzale fax: +386 (0)1 72 17 888
Slovenia, Europe
--
One must learn by doing the thing; for though you think you know it,
you have no certainty until you try. Sophocles ~ 450 B.C.
--
-Original Message-
From: Duncan Murdoch [mailto:[EMAIL PROTECTED]
Sent: sre 2005-06-08 15:53
To: Gorjanc Gregor
Cc: r-help@stat.math.ethz.ch
Subject: Re: [R] Random seed problem in MCMC coupling of chains
On 6/8/2005 9:27 AM, Gorjanc Gregor wrote:
Hello!
I am performing coupling of chains in MCMC and I need the same value
of seed for two chains. I will show demo of what I want:
R code, which might show my example is:
niter - 3
nchain - 2
tmpSeed - 123
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(tmpSeed)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
tmpSeed - .Random.seed
}
}
I get this:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.43588
iter: 2 chain: 2 runif: 0.43588
iter: 3 chain: 1 runif: 0.43588
iter: 3 chain: 2 runif: 0.43588
but I would like to get:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.67676
iter: 2 chain: 2 runif: 0.67676
iter: 3 chain: 1 runif: 0.12368
iter: 3 chain: 2 runif: 0.12368
Note that seed value is of course changing, but it is parallel
between chains.
set.seed takes an integer, but .Random.seed is a complicated vector.
You need to play with .Random.seed directly, and move your setting of
tmpSeed out of the inner loop, i.e.
niter - 3
nchain - 2
set.seed(123)
tmpSeed - .Random.seed
for (i in 1:niter) { # iterations
+ for (j in 1:nchain) { # chains
+ .Random.seed - tmpSeed
+ a - runif(1)
+ cat(iter:, i, chain:, j, runif:, a, \n)
+ }
+ tmpSeed - .Random.seed
+ }
iter: 1 chain: 1 runif: 0.2875775
iter: 1 chain: 2 runif: 0.2875775
iter: 2 chain: 1 runif: 0.7883051
iter: 2 chain: 2 runif: 0.7883051
iter: 3 chain: 1 runif: 0.4089769
iter: 3 chain: 2 runif: 0.4089769
However, heed the warnings in ?set.seed: with some generators
.Random.seed *does not* contain the full state of the RNG. As far as I
know there is no way to obtain the full state.
Duncan Murdoch
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Solve f(x) = 0
Hi! I´m need a function that solves the equation f(x) = 0 (i.e. the root of the function) when f is a nonlinear function. Is there any? I´ve tried nlm and optim on the square of the function but the solution is very unstable. Thanks before hand. / Fredrik Thuring -- This e-mail and any attachment may be confidential and may also be privileged. If you are not the intended recipient, please notify us immediately and then delete this e-mail and any attachment without retaining copies or disclosing the contents thereof to any other person. Thank you. -- [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] variable png-size in multiple figure environment
On Wed, 2005-06-08 at 16:46 +0200, Claude Thato Kobe wrote: If I put only one plot in a png, I can drag the generated png on the edges and it changes its size by this. But when I put several plots on the same sheet using mfrow, the size can no longer be changed when viewing the file, and the resolution is bad. What do I need to do to keep the variability of size in a multiple figure environment? I have no problem resizing either type of PNG file using ImageMagick under FC3 by dragging a corner with a mouse, so I suspect the behavior may be a function of your viewing application. For example, using other viewers I can zoom in and out, but not drag to resize. Don't use a bitmapped graphic format if you need to resize, as these do not resize well, especially diagonal and curved lines. You should use a vector based format if you need to re-size the image and maintain image quality. Depending upon your platform (which you do not indicate) use: WMF PS/EPS PDF WMF will be available under Windows. The choice may depend upon what you intend to do with the plots. Alternatively, pre-specify the graphic device height and width using the appropriate arguments to the functions so that the output image is as close as possible to the actual size that you require. Note that for bitmapped plots, these will be in pixels for the png() function and inches if you use the bitmap() function, which requires ghostscript. The combination of pixels and resolution settings will determine linear size for a given output device when using bitmapped graphics. See the help for the specific devices you wish to use or ?Devices for general help. HTH, Marc Schwartz __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] get level combinations from by list
I think these are not stored by by() and they are produced inside
print.by(). One way to get them could be:
dn - dimnames(grouped)
dn - expand.grid(dn)
dn - dn[!sapply(grouped[1:nrow(dn)], is.null), ]
dn
rownames(dfr) - apply(dn, 1, paste, collapse = -)
dfr
I hope it helps.
Best,
Dimitris
Dimitris Rizopoulos
Ph.D. Student
Biostatistical Centre
School of Public Health
Catholic University of Leuven
Address: Kapucijnenvoer 35, Leuven, Belgium
Tel: +32/16/336899
Fax: +32/16/337015
Web: http://www.med.kuleuven.ac.be/biostat/
http://www.student.kuleuven.ac.be/~m0390867/dimitris.htm
- Original Message -
From: bogdan romocea [EMAIL PROTECTED]
To: R-help@stat.math.ethz.ch
Sent: Wednesday, June 08, 2005 4:56 PM
Subject: [R] get level combinations from by list
Dear useRs,
Given this code I end up with a list of class by:
a - sample(1:5,200,replace=TRUE)
b - sample(c(v1,v2,v3),200,replace=TRUE)
c - sample(c(11,22,33),200,replace=TRUE)
data - runif(200)
grouped - by(data,list(a,b,c),function(x) {c(min=min(x),max=max(x),
median=round(median(x),digits=2),mean=round(mean(x),digits=2))})
dfr - do.call(rbind,grouped)#the levels are missing
#--
grouped
typeof(grouped)
class(grouped)
dimnames(grouped)
How do I get at the levels of the 'group by' variables for each
subset? For example, from this part of the by list I want 4, v2
and
33:
: 4
: v2
: 33
min maxmedian mean
0.3897450 0.9215315 0.730 0.670
---
Thank you,
b.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Solve f(x) = 0
On Wed, 8 Jun 2005 17:08:35 +0200 Fredrik Thuring wrote: Hi! I´m need a function that solves the equation f(x) = 0 (i.e. the root of the function) when f is a nonlinear function. Is there any? I´ve tried nlm and optim on the square of the function but the solution is very unstable. Look at ?uniroot hth, Z Thanks before hand. / Fredrik Thuring - - This e-mail and any attachment may be confidential and may also be privileged. If you are not the intended recipient, please notify us immediately and then delete this e-mail and any attachment without retaining copies or disclosing the contents thereof to any other person. Thank you. - - [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Solve f(x) = 0
On Wed, 2005-06-08 at 17:08 +0200, Fredrik Thuring wrote: Hi! Im need a function that solves the equation f(x) = 0 (i.e. the root of the function) when f is a nonlinear function. Is there any? Ive tried nlm and optim on the square of the function but the solution is very unstable. Your colleague Martin Englund asked the same question. If x is scalar, you can use uniroot(). If it's a polynomial then polyroot() finds all the zeros. I hope that helps. Martyn __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Solve f(x) = 0
Have you plotted f(x)? Plot(s) might help reveal why the solution is very unstable. If f is a function of a univariate x, this is trivial. If f is a function of a bivariate x, use something like contour or persp. If higher dimensions, I might use something like fit - optim(... hessian=TRUE) eigen(fit$hessian) hope this helps. spencer graves Martyn Plummer wrote: On Wed, 2005-06-08 at 17:08 +0200, Fredrik Thuring wrote: Hi! Im need a function that solves the equation f(x) = 0 (i.e. the root of the function) when f is a nonlinear function. Is there any? Ive tried nlm and optim on the square of the function but the solution is very unstable. Your colleague Martin Englund asked the same question. If x is scalar, you can use uniroot(). If it's a polynomial then polyroot() finds all the zeros. I hope that helps. Martyn __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Random seed problem in MCMC coupling of chains
That could be addressed like this (where changing the offset
changes the experiment).
offset - 123
niter - 3
nchain - 2
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(i+offset)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
}
}
On 6/8/05, Paul Gilbert [EMAIL PROTECTED] wrote:
Beware that your easy trick will give you the same result every time you
run it. You need a better scheme if you actually intend to get a new
experiment each time you run it.
Paul
Gorjanc Gregor wrote:
Thanks to Duncan, Dimitris as well as James for answers. I'll provide
here also example from James, which seems to be the easiest of them
all and was not posted to the list:
niter - 3
nchain - 2
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(i)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
}
}
Note that seed is set with iteration counter. This is really neat and
simple. I am just wondering if this is OK from RNG point of view. Can
someone comment on that?
Lep pozdrav / With regards,
Gregor Gorjanc
--
University of Ljubljana
Biotechnical FacultyURI: http://www.bfro.uni-lj.si/MR/ggorjan
Zootechnical Department mail: gregor.gorjanc at bfro.uni-lj.si
Groblje 3 tel: +386 (0)1 72 17 861
SI-1230 Domzale fax: +386 (0)1 72 17 888
Slovenia, Europe
--
One must learn by doing the thing; for though you think you know it,
you have no certainty until you try. Sophocles ~ 450 B.C.
--
-Original Message-
From: Duncan Murdoch [mailto:[EMAIL PROTECTED]
Sent: sre 2005-06-08 15:53
To: Gorjanc Gregor
Cc: r-help@stat.math.ethz.ch
Subject: Re: [R] Random seed problem in MCMC coupling of chains
On 6/8/2005 9:27 AM, Gorjanc Gregor wrote:
Hello!
I am performing coupling of chains in MCMC and I need the same value
of seed for two chains. I will show demo of what I want:
R code, which might show my example is:
niter - 3
nchain - 2
tmpSeed - 123
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(tmpSeed)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
tmpSeed - .Random.seed
}
}
I get this:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.43588
iter: 2 chain: 2 runif: 0.43588
iter: 3 chain: 1 runif: 0.43588
iter: 3 chain: 2 runif: 0.43588
but I would like to get:
iter: 1 chain: 1 runif: 0.43588
iter: 1 chain: 2 runif: 0.43588
iter: 2 chain: 1 runif: 0.67676
iter: 2 chain: 2 runif: 0.67676
iter: 3 chain: 1 runif: 0.12368
iter: 3 chain: 2 runif: 0.12368
Note that seed value is of course changing, but it is parallel
between chains.
set.seed takes an integer, but .Random.seed is a complicated vector.
You need to play with .Random.seed directly, and move your setting of
tmpSeed out of the inner loop, i.e.
niter - 3
nchain - 2
set.seed(123)
tmpSeed - .Random.seed
for (i in 1:niter) { # iterations
+ for (j in 1:nchain) { # chains
+ .Random.seed - tmpSeed
+ a - runif(1)
+ cat(iter:, i, chain:, j, runif:, a, \n)
+ }
+ tmpSeed - .Random.seed
+ }
iter: 1 chain: 1 runif: 0.2875775
iter: 1 chain: 2 runif: 0.2875775
iter: 2 chain: 1 runif: 0.7883051
iter: 2 chain: 2 runif: 0.7883051
iter: 3 chain: 1 runif: 0.4089769
iter: 3 chain: 2 runif: 0.4089769
However, heed the warnings in ?set.seed: with some generators
.Random.seed *does not* contain the full state of the RNG. As far as I
know there is no way to obtain the full state.
Duncan Murdoch
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Specifying medoids in PAM?
David == David Finlayson [EMAIL PROTECTED]
on Wed, 8 Jun 2005 09:24:54 -0700 writes:
David Sorry, I wasn't trying to submit a bug report just yet.
the posting guide asks you to provide reproducible examples, in
any case, not just for bug reports ...
{and strictly speaking, you still haven't provided one, since
it's a bit painful to read in your table below -- because of the
extra row names ... but here I'm nit picking a bit }
David I wanted to see if I was using the command correctly.
Yes, you were.
David Here is a complete example (the row names have numbers in them):
library(cluster)
stats.table
David slope.for slope.ter prof.len break.elv
David 1 ALLYN W 0.09 0.05 63.3 1.46
David 2 ARCADIA 0.12 0.09 40.8 0.43
David 3 BURFOOT 0.12 0.03 58.9 0.78
David 4 CAMA BE 0.15 0.05 50.4 0.88
David 5 DASH PO 0.09 0.01290.5 3.14
David 6 EAGLE C 0.09 0.01226.4 1.49
..
pam(stats.table, metric=euclidean, stand=TRUE, medoids=c(1,3,20,2,5),
k=5)
David This command crashes RGUI.exe and windows sends an error report to
David Microsoft. It also crashes if I first subtract the NA rows from
David stats.table.
I can confirm to get segmentation faults using this example data
with k=5 , so effectively, it seems you've uncovered a bug in pam().
I will investigate and patch eventually.
Thank you for the report.
Martin Maechler,
ETH Zurich
David Thanks,
David David
David On 6/8/05, Martin Maechler [EMAIL PROTECTED] wrote:
David == David Finlayson [EMAIL PROTECTED]
on Tue, 7 Jun 2005 12:11:25 -0700 writes:
David I am using the PAM algorithm in the CLUSTER library.
David When I allow PAM to seed the medoids using the default __build__
David algorithm things work
David well:
pam(stats.table, metric=euclidean, stand=TRUE, k=5)
David But I have some clusters from a Hierarchical analysis that I would
David like to use as seeds for the PAM algorithm. I can't figure what the
David mediod argument wants. When I put in the five integer indicies for
the
David observations in stats.table that I would like to use as seeds (the
row
David numbers), I segfault R.
pam(stats.table, metric=euclidean, stand=TRUE,
medoids=c(1,3,20,2,5), k=5)
David *** R Crashes ***
this is not very helpful.
Can you please READ the posting guide and then
do as the guide says :
- post a *reproducible* example
- tell more about what happens when R Crashes
David Here is my version info:
version
David _
David platform i386-pc-mingw32
David arch i386
David os mingw32
David system i386, mingw32
David status
David major2
David minor0.1
David year 2004
David month11
David day 15
David language R
David Any guidance would be appreciated.
David David
David --
David David Finlayson
David Marine Geology Geophysics
David School of Oceanography
David Box 357940
David University of Washington
David Seattle, WA 98195-7940
David USA
David Office: Marine Sciences Building, Room 112
David Phone: (206) 616-9407
David Web: http://students.washington.edu/dfinlays
David __
David R-help@stat.math.ethz.ch mailing list
David https://stat.ethz.ch/mailman/listinfo/r-help
David PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
David --
David David Finlayson
David Marine Geology Geophysics
David School of Oceanography
David Box 357940
David University of Washington
David Seattle, WA 98195-7940
David USA
David Office: Marine Sciences Building, Room 112
David Phone: (206) 616-9407
David Web: http://students.washington.edu/dfinlays
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] Random seed problem in MCMC coupling of chains
Thanks to Paul and Gabor for additional tips/examples. Actually, I find
Pauls suggestion with setRNG also nice and is exactly what I wanted.
Paul, if I understand this correctly, your suggestion with setRNG does not
alter RNG flow, it just takes care that chains really have equal seeds.
I remember that I have read somewhere that destroying RNG flow over and
over to get real randomness is not a good idea. Can someone confirm this?
niter - 3
nchain - 2
for (i in 1:niter) { # iterations
tmpSeed - setRNG()
for (j in 1:nchain) { # chains
setRNG(tmpSeed)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
}
}
iter: 1 chain: 1 runif: 0.8160078
iter: 1 chain: 2 runif: 0.8160078
iter: 2 chain: 1 runif: 0.4909793
iter: 2 chain: 2 runif: 0.4909793
iter: 3 chain: 1 runif: 0.4425924
iter: 3 chain: 2 runif: 0.4425924
[... removed other stuff ...]
Lep pozdrav / With regards,
Gregor Gorjanc
--
University of Ljubljana
Biotechnical FacultyURI: http://www.bfro.uni-lj.si/MR/ggorjan
Zootechnical Department mail: gregor.gorjanc at bfro.uni-lj.si
Groblje 3 tel: +386 (0)1 72 17 861
SI-1230 Domzale fax: +386 (0)1 72 17 888
Slovenia, Europe
--
One must learn by doing the thing; for though you think you know it,
you have no certainty until you try. Sophocles ~ 450 B.C.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to run a script at the beginning of an interacive session ?
On 6/8/05, Marc Schwartz [EMAIL PROTECTED] wrote:
snips
Thanks Marc,
Cool ideas, thanks! Building on them, here is a twist.
With this line in .Rprofile:
.First - function() {if(Sys.getenv(R_PROJECT)!=)
source(Sys.getenv(R_PROJECT))}
R can be run with this bash command line:
## export R_PROJECT=project_A.R; R; unset R_PROJECT
Alternatively, R could be run with this bash command line:
## wrapR project_4.R
contents of wrapR
#! /bin/bash
export R_PROJECT=$1
R
end
What do you think?
Thanks in advance,
Mike
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Random seed problem in MCMC coupling of chains
Here is a small variation. We define a list to hold
the last seed for each chain. Each time we enter the simulation
for a chain we use that seed and each time we exit we update it.
The loop becomes simpler since the setup is all done prior
to looping and everything else is done in the inner loop.
Note that a double loop with nothing between the first and second
for is really like a single loop over the i,j pairs so its presumably
easier to understand.
library(setRNG)
set.seed(123)
niter - 3; nchain - 2
chain - lapply(1:nchain, function(x) setRNG())
for(i in 1:niter)
for(j in 1:nchain) {
setRNG(chain[[j]]) # get seed
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
chain[[j]] - setRNG() # save seed
}
On 6/8/05, Gorjanc Gregor [EMAIL PROTECTED] wrote:
Thanks to Paul and Gabor for additional tips/examples. Actually, I find
Pauls suggestion with setRNG also nice and is exactly what I wanted.
Paul, if I understand this correctly, your suggestion with setRNG does not
alter RNG flow, it just takes care that chains really have equal seeds.
I remember that I have read somewhere that destroying RNG flow over and
over to get real randomness is not a good idea. Can someone confirm this?
niter - 3
nchain - 2
for (i in 1:niter) { # iterations
tmpSeed - setRNG()
for (j in 1:nchain) { # chains
setRNG(tmpSeed)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
}
}
iter: 1 chain: 1 runif: 0.8160078
iter: 1 chain: 2 runif: 0.8160078
iter: 2 chain: 1 runif: 0.4909793
iter: 2 chain: 2 runif: 0.4909793
iter: 3 chain: 1 runif: 0.4425924
iter: 3 chain: 2 runif: 0.4425924
[... removed other stuff ...]
Lep pozdrav / With regards,
Gregor Gorjanc
--
University of Ljubljana
Biotechnical FacultyURI: http://www.bfro.uni-lj.si/MR/ggorjan
Zootechnical Department mail: gregor.gorjanc at bfro.uni-lj.si
Groblje 3 tel: +386 (0)1 72 17 861
SI-1230 Domzale fax: +386 (0)1 72 17 888
Slovenia, Europe
--
One must learn by doing the thing; for though you think you know it,
you have no certainty until you try. Sophocles ~ 450 B.C.
--
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Random seed problem in MCMC coupling of chains
On Wed, Jun 08, 2005 at 12:55:07PM -0400, Gabor Grothendieck wrote:
That could be addressed like this (where changing the offset
changes the experiment).
offset - 123
niter - 3
nchain - 2
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(i+offset)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
}
}
On 6/8/05, Paul Gilbert [EMAIL PROTECTED] wrote:
Beware that your easy trick will give you the same result every time you
run it. You need a better scheme if you actually intend to get a new
experiment each time you run it.
That's not a bad thing per se; in fact, it's a good thing to be able
to exactly reproduce the results you obtained with your software.
Personally, I dislike the convenience feature of random number
generators of generating a seed, frequently based on the system time,
if none has been set explicitly; I always set a seed and frequently
make the seed a commandline option or part of the control parameter
file or the like.
From this perspective, Gabor's solution seems perfect to me.
Paul
Gorjanc Gregor wrote:
Thanks to Duncan, Dimitris as well as James for answers. I'll provide
here also example from James, which seems to be the easiest of them
all and was not posted to the list:
niter - 3
nchain - 2
for (i in 1:niter) { # iterations
for (j in 1:nchain) { # chains
set.seed(i)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
}
}
Note that seed is set with iteration counter. This is really neat and
simple. I am just wondering if this is OK from RNG point of view. Can
someone comment on that?
The only concern I could think about is the case of a bad random number
generator, in which the first couple of values are not entirely
uncorrelated to the seed. But I'd be very surprised if that was a
problem with R's RNGs -- I guess it's memories of lousy implementations
C library rand() functions that make me write this remark.
Best regards, Jan
--
+- Jan T. Kim ---+
|*NEW*email: [EMAIL PROTECTED] |
|*NEW*WWW: http://www.cmp.uea.ac.uk/people/jtk |
*-= hierarchical systems are for files, not for humans =-*
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] install package snow
Hi, I have difficulty in installing package snow from CRAN. Somehow, this package is not shown up in the available package list when I tried to use the GUI interface. So I have to go to CRAN to download the source file: snow_0.2-1.tar.gz. Then I tried to install it using the following command: install.packages(repos=NULL, pkgs='C:\Documents and Settings\Desktop\snow_0.2-1.tar.gz', type='source', lib = 'c:/program files/r/rw2010/library/') Warning message: installation of package 'C:Documents and Settingsyzhang24Desktopsnow_0.2-1.tar.gz' had non-zero exit status in: install.packages(repos = NULL, pkgs = C:Documents and Settingsyzhang24Desktopsnow_0.2-1.tar.gz, What's the problem here, and how should I handle it? Thanks Steve __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Course***R/S-plus Programming at 4 locations - July 2005
XLSolutions Corporation (www.xlsolutions-corp.com) is proud to announce 2-day R/S-plus Fundamentals and Programming Techniques at 4 locations nationwide. www.xlsolutions-corp.com/courselist.htm Seattle, WA - July 11th-12th, 2005 San Francisco -- July 14th-15th, 2005 New York, NY -- July 25th - 26th,2005 Boston, MA - July 28th - 29th, 2005 Reserve your seat now at the early bird rates! Payment due AFTER the class Course Description: This two-day beginner to intermediate R/S-plus course focuses on a broad spectrum of topics, from reading raw data to a comparison of R and S. We will learn the essentials of data manipulation, graphical visualization and R/S-plus programming. We will explore statistical data analysis tools,including graphics with data sets. How to enhance your plots, build your own packages (librairies) and connect via ODBC,etc. We will perform some statistical modeling and fit linear regression models. Participants are encouraged to bring data for interactive sessions With the following outline: - An Overview of R and S - Data Manipulation and Graphics - Using Lattice Graphics - A Comparison of R and S-Plus - How can R Complement SAS? - Writing Functions - Avoiding Loops - Vectorization - Statistical Modeling - Project Management - Techniques for Effective use of R and S - Enhancing Plots - Using High-level Plotting Functions - Building and Distributing Packages (libraries) - Connecting; ODBC, Rweb, Orca via sockets and via Rjava Interested in R/Splus Advanced course? email us. Email us for group discounts. Email Sue Turner: [EMAIL PROTECTED] Phone: 206-686-1578 Visit us: www.xlsolutions-corp.com/training.htm Please let us know if you and your colleagues are interested in this classto take advantage of group discount. Register now to secure your seat! Interested in R/Splus Advanced course? email us. Cheers, Elvis Miller, PhD Manager Training. XLSolutions Corporation 206 686 1578 www.xlsolutions-corp.com [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] install package snow
I think for 'snow' you need to download the source and build it yourself. I assume you are on Windows (you didn't say!). 'snow' depends on external programs so it would seem unwise to have a binary package on CRAN. -roger Steve Adams wrote: Hi, I have difficulty in installing package snow from CRAN. Somehow, this package is not shown up in the available package list when I tried to use the GUI interface. So I have to go to CRAN to download the source file: snow_0.2-1.tar.gz. Then I tried to install it using the following command: install.packages(repos=NULL, pkgs='C:\Documents and Settings\Desktop\snow_0.2-1.tar.gz', type='source', lib = 'c:/program files/r/rw2010/library/') Warning message: installation of package 'C:Documents and Settingsyzhang24Desktopsnow_0.2-1.tar.gz' had non-zero exit status in: install.packages(repos = NULL, pkgs = C:Documents and Settingsyzhang24Desktopsnow_0.2-1.tar.gz, What's the problem here, and how should I handle it? Thanks Steve __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Roger D. Peng http://www.biostat.jhsph.edu/~rpeng/ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] install package snow
Please do read the FAQ Q7.8. Once you have followed that, it should install from the sources (it does for me). However, to make it work you need other things which are non-trivial to get working on Windows: see http://www.stat.uiowa.edu/~luke/R/cluster/cluster.html for the lower-level interfaces. I think snow over PVM has been reported to work under Windows, but I don't know about the others. There is no binary version of the package, as testing that would be tedious at best. On Wed, 8 Jun 2005, Steve Adams wrote: Hi, I have difficulty in installing package snow from CRAN. Somehow, this package is not shown up in the available package list when I tried to use the GUI interface. So I have to go to CRAN to download the source file: snow_0.2-1.tar.gz. Then I tried to install it using the following command: install.packages(repos=NULL, pkgs='C:\Documents and Settings\Desktop\snow_0.2-1.tar.gz', type='source', lib = 'c:/program files/r/rw2010/library/') Warning message: installation of package 'C:Documents and Settingsyzhang24Desktopsnow_0.2-1.tar.gz' had non-zero exit status in: install.packages(repos = NULL, pkgs = C:Documents and Settingsyzhang24Desktopsnow_0.2-1.tar.gz, What's the problem here, and how should I handle it? Thanks Steve __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Brian D. Ripley, [EMAIL PROTECTED] 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@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] install package snow
Le 08.06.2005 20:51, Steve Adams a écrit : Hi, I have difficulty in installing package snow from CRAN. Somehow, this package is not shown up in the available package list when I tried to use the GUI interface. So I have to go to CRAN to download the source file: snow_0.2-1.tar.gz. Then I tried to install it using the following command: install.packages(repos=NULL, pkgs='C:\Documents and Settings\Desktop\snow_0.2-1.tar.gz', type='source', lib = 'c:/program files/r/rw2010/library/') Double backslash are required here, as in : install.packages(repos=NULL, pkgs='C:\\Documents and Settings\\Desktop\\snow_0.2-1.tar.gz', type='source',lib = 'c:/program files/r/rw2010/library/') (but the tools are required if you want to build from source on windows, so that may not solve the problem) Romain Warning message: installation of package 'C:Documents and Settingsyzhang24Desktopsnow_0.2-1.tar.gz' had non-zero exit status in: install.packages(repos = NULL, pkgs = C:Documents and Settingsyzhang24Desktopsnow_0.2-1.tar.gz, What's the problem here, and how should I handle it? Thanks Steve -- visit the R Graph Gallery : http://addictedtor.free.fr/graphiques ~ ~~ Romain FRANCOIS - http://addictedtor.free.fr ~~ Etudiant ISUP - CS3 - Industrie et Services ~~http://www.isup.cicrp.jussieu.fr/ ~~ Stagiaire INRIA Futurs - Equipe SELECT ~~ http://www.inria.fr/recherche/equipes/select.fr.html~~ ~ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
FW: [R] Random seed problem in MCMC coupling of chains
And a last post from Paul Gilbert. Thanks to all! This disscusion was
really beneficial for me!
-Original Message-
From: Paul Gilbert [mailto:[EMAIL PROTECTED]
Sent: sre 2005-06-08 21:01
To: Gorjanc Gregor
Subject: Re: [R] Random seed problem in MCMC coupling of chains
Gorjanc Gregor wrote:
Thanks to Paul and Gabor for additional tips/examples. Actually, I find
Pauls suggestion with setRNG also nice and is exactly what I wanted.
Paul, if I understand this correctly, your suggestion with setRNG does not
alter RNG flow, it just takes care that chains really have equal seeds.
I remember that I have read somewhere that destroying RNG flow over and
over to get real randomness is not a good idea. Can someone confirm this?
In general it is a bad idea to make up your own scheme for setting or
resetting the RNG. People put a lot of work into studying the properties
of a RNG. When you mess with it then it is unclear what the result will
be. It certainly won't be tested unless you test it yourself. If your
intention is to do research on RNGs then you may want to do that, but if
your intention is to do other research and just use the RNG, then don't
mess with it by resetting it with your own scheme.
One additional thing you may want to do is record the initial setting of
the RNG information so that you can reproduce the experiment if you want
to (see modification below). The idea in setRNG is to not interfere with
the flow, only add a few utilities to help record and reset everything
when that is what is required.
In your example it is important that you generate the same number of
random numbers in each pass through the chain. If that is not the case
then even with the setRNG utilities there is a subtle change that you
are introducing.
HTH,
Paul
niter - 3
nchain - 2
startingRNG - setRNG()
for (i in 1:niter) { # iterations
tmpSeed - setRNG()
for (j in 1:nchain) { # chains
setRNG(tmpSeed)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
}
}
iter: 1 chain: 1 runif: 0.8160078
iter: 1 chain: 2 runif: 0.8160078
iter: 2 chain: 1 runif: 0.4909793
iter: 2 chain: 2 runif: 0.4909793
iter: 3 chain: 1 runif: 0.4425924
iter: 3 chain: 2 runif: 0.4425924
[... removed other stuff ...]
Lep pozdrav / With regards,
Gregor Gorjanc
--
University of Ljubljana
Biotechnical FacultyURI: http://www.bfro.uni-lj.si/MR/ggorjan
Zootechnical Department mail: gregor.gorjanc at bfro.uni-lj.si
Groblje 3 tel: +386 (0)1 72 17 861
SI-1230 Domzale fax: +386 (0)1 72 17 888
Slovenia, Europe
--
One must learn by doing the thing; for though you think you know it,
you have no certainty until you try. Sophocles ~ 450 B.C.
--
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Random seed problem in MCMC coupling of chains
Gorjanc Gregor wrote:
Thanks to Paul and Gabor for additional tips/examples. Actually, I find
Pauls suggestion with setRNG also nice and is exactly what I wanted.
Paul, if I understand this correctly, your suggestion with setRNG does not
alter RNG flow, it just takes care that chains really have equal seeds.
I remember that I have read somewhere that destroying RNG flow over and
over to get real randomness is not a good idea. Can someone confirm this?
That's in Brian Ripley's Simulation book, and certainly in other places.
Kjetil
niter - 3
nchain - 2
for (i in 1:niter) { # iterations
tmpSeed - setRNG()
for (j in 1:nchain) { # chains
setRNG(tmpSeed)
a - runif(1)
cat(iter:, i, chain:, j, runif:, a, \n)
}
}
iter: 1 chain: 1 runif: 0.8160078
iter: 1 chain: 2 runif: 0.8160078
iter: 2 chain: 1 runif: 0.4909793
iter: 2 chain: 2 runif: 0.4909793
iter: 3 chain: 1 runif: 0.4425924
iter: 3 chain: 2 runif: 0.4425924
[... removed other stuff ...]
Lep pozdrav / With regards,
Gregor Gorjanc
--
University of Ljubljana
Biotechnical FacultyURI: http://www.bfro.uni-lj.si/MR/ggorjan
Zootechnical Department mail: gregor.gorjanc at bfro.uni-lj.si
Groblje 3 tel: +386 (0)1 72 17 861
SI-1230 Domzale fax: +386 (0)1 72 17 888
Slovenia, Europe
--
One must learn by doing the thing; for though you think you know it,
you have no certainty until you try. Sophocles ~ 450 B.C.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
--
Kjetil Halvorsen.
Peace is the most effective weapon of mass construction.
-- Mahdi Elmandjra
--
No virus found in this outgoing message.
Checked by AVG Anti-Virus.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: FW: [R] Random seed problem in MCMC coupling of chains
On 6/8/05, Gorjanc Gregor [EMAIL PROTECTED] wrote: And a last post from Paul Gilbert. -Original Message- From: Paul Gilbert [mailto:[EMAIL PROTECTED] In your example it is important that you generate the same number of random numbers in each pass through the chain. If that is not the case then even with the setRNG utilities there is a subtle change that you are introducing. Note that this is actually one of the advantages of the last solution I posted, namely the one with the chain list. It maintains the flow along each chain even if different chains use different number of calls to the random number generator. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to run a script at the beginning of an interacive session ?
On Wed, 2005-06-08 at 11:32 -0700, Mike R wrote:
On 6/8/05, Marc Schwartz [EMAIL PROTECTED] wrote:
snips
Thanks Marc,
Cool ideas, thanks! Building on them, here is a twist.
With this line in .Rprofile:
.First - function() {if(Sys.getenv(R_PROJECT)!=)
source(Sys.getenv(R_PROJECT))}
R can be run with this bash command line:
## export R_PROJECT=project_A.R; R; unset R_PROJECT
Alternatively, R could be run with this bash command line:
## wrapR project_4.R
contents of wrapR
#! /bin/bash
export R_PROJECT=$1
R
end
What do you think?
Thanks in advance,
Mike
Mike,
At this point, I think whatever approach to implementing and managing
the process you might find easier is the way to go. Whether it be the
above or my approach.
The initial iteration resulted in overwriting the default ~/.Rprofile
file, which could cause problems if you had other options/settings in
it.
The key is ease of implementation, reduction in conflicts/errors and
ongoing maintenance if you find yourself needing multiple startup
options.
HTH,
Marc
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] [R-pkgs] New CRAN package dyn
dyn is an R package that facilitates the use of regression using time series data with lags and diffs (known as dynamic regression). It is a lightweight package that has no facilities of its own but leverages off the various time series and regression functions in R to make it easier to use them together. Its features include: - many regression functions. It can be used with lm, glm, loess, rq, randomForest, lqs, rlm and any other regression functions that use model.frame and are written in the style of lm. - many time series classes. It can be used with ts, zooreg, zoo, its, and irts time series classes. This covers regular, weakly regular and irregular time series classes. - missing values. Time series may have missing values including internal missing values. Both na.omit and na.exclude are supported. - good citizen. It does not replace the regression functions but rather works with them by providing new methods to the standard R generics: model.frame, resid, fitted, predict, update, anova and $. - ease of use. dyn enables one to use the same regression functions (lm, glm, etc.) using the same syntax one has always used. Just preface the regression function name with dyn$ and it is transformed into a regression function that can handle time series: dyn$lm( y ~ x + lag(x) + diff(w) ) # lm dyn$loess( y ~ x + lag(x) + diff(w) ) # loess - modular. dyn can be used with any regression function that uses model.frame and is written in the style of lm. Additional classes can be added to dyn simply by adding new methods. dyn is modular so such updates can be made without changing dyn, itself. - documentation. It includes a help page and six demos. ?dyn # help file demo() # look under dyn for list of demos demo(dyn-rq) # runs indicated dyn demo The package is available on CRAN. Comments/questions welcome. ___ R-packages mailing list [EMAIL PROTECTED] https://stat.ethz.ch/mailman/listinfo/r-packages __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Specifying medoids in PAM?
MM == Martin Maechler [EMAIL PROTECTED]
on Wed, 8 Jun 2005 18:57:55 +0200 writes:
David == David Finlayson [EMAIL PROTECTED]
on Wed, 8 Jun 2005 09:24:54 -0700 writes:
David Sorry, I wasn't trying to submit a bug report just yet.
MM the posting guide asks you to provide reproducible examples, in
MM any case, not just for bug reports ...
MM {and strictly speaking, you still haven't provided one, since
MM it's a bit painful to read in your table below -- because of the
MM extra row names ... but here I'm nit picking a bit }
David I wanted to see if I was using the command correctly.
MM Yes, you were.
pam(stats.table, metric=euclidean, stand=TRUE, medoids=c(1,3,20,2,5),
k=5)
David This command crashes RGUI.exe and windows sends an error report to
David Microsoft. It also crashes if I first subtract the NA rows from
David stats.table.
MM I can confirm to get segmentation faults using this example data
MM with k=5 , so effectively, it seems you've uncovered a bug in pam().
MM I will investigate and patch eventually.
I found and fixed the bug:
Some part of the C code was assuming that the indices in
'medoids' were sorted (increasingly).
I.e., for the moment you can easily work around the problem by
using
pam(stats.table, , medoids=c(1,2,3,5,20), k=5)
instead of
pam(stats.table, , medoids=c(1,3,20,2,5), k=5)
The next version of the cluster package which allows to specify
the fuzzyness exponent in fanny() will have this problem
fixed.
Martin Maechler,
ETH Zurich
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Feedback requested on a solution to specifying R-scripts on the command line
Hi All,
I'm new to R and would like to get started on the right
foot, so to say. So I am requesting feedback on a
paradigm for organizing R projects. In esence, I am
trying to organize my R projects in the same way I
organize my numerical simulations. But it feels like
I am imposing
In particular, I'm not sure why the ability to specify
R-scripts on the command line is not built in to R.
For example, am I not understanding the traditional
paradigm for R projects? If so, is there documentation
that gives a practical orientation to R projects?
From the command line, I'd like to be able to start an
interactive R session and at the same time, be able
to specify on the command line a project-specific or
task-specific script (R code) that is to be executed
at the beginning of the session.
My current solution is the following:
Add the following line to dot-Rprofile:
.First - function()
{
if( Sys.getenv(R_CUSTOM) != )
{
for ( filename in strsplit(Sys.getenv(R_CUSTOM), )[[1]] )
{
source( filename )
}
}
}
Then set the environement variable R_CUSTOM using a wrapper.
Here is a wrapper for bash:
#! /bin/bash
export R_CUSTOM=$*
R
If the wrapper file is named wrapR then here are some sample
command lines (my shell prompt is set to ##):
## wrapR
## wrapR project_kit.R
## wrapR project_kit.R 00task_verify_mikes_datafile.R
## wrapR project_kit.R 02task_verify_jeffs_datafile.R
## wrapR project_kit.R 11task_analyze_pooled_data.R
## wrapR project_kit.R 21task_generate_figures_for_executivesummary.R
## wrapR project_kit.R 22task_generate_figures_for_paper.R
Normally these could be run using redirection,
as in the following command line:
## R 22kit_n_task_generate_figures_for_paper.R
but in that case, R exits when through.
Any comments would be welcome.
Thanks in advance,
Mike
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] How to plot more than 3 sets in Venn Diagrams?
I'm trying to plot Venn diagrams with more than 3 sets (5 actually) in order to describe graphically the genetic variation between populations. I tried the limma library but realised it can only plot 3 sets. Is there any solution? Of course I could plot the chart manually but it'll take too long (have other datasets). One of my dataset is given below. THanks for any advice. j AA CH EA IN MY [1,] 0 0 0 1 0 [2,] 1 0 0 0 0 [3,] 1 0 0 0 0 [4,] 0 0 0 0 1 [5,] 1 0 0 0 0 [6,] 1 0 0 0 0 [7,] 1 0 0 0 0 [8,] 1 0 0 1 0 [9,] 1 0 0 0 0 [10,]1 0 0 0 0 [11,]1 0 0 0 0 [12,]0 0 0 0 1 [13,]1 0 0 0 0 [14,]0 0 1 0 0 [15,]1 0 0 0 0 [16,]0 0 1 0 0 [17,]1 0 0 0 0 [18,]0 0 0 1 0 [19,]1 0 0 1 0 [20,]0 0 1 0 0 [21,]0 1 0 0 0 [22,]1 0 0 0 0 [23,]0 0 1 0 0 [24,]1 0 1 1 1 [25,]0 1 0 0 0 [26,]1 0 0 0 0 [27,]1 0 0 0 0 [28,]0 0 0 1 0 [29,]0 0 0 1 0 [30,]1 0 0 0 0 [31,]1 0 0 0 0 [32,]0 0 0 1 0 [33,]0 0 0 1 0 [34,]0 1 0 0 0 [35,]1 0 0 0 0 [36,]0 0 0 1 0 [37,]0 0 1 1 0 [38,]1 0 0 1 0 [39,]0 0 0 1 0 [40,]0 0 0 1 0 [41,]0 1 0 0 0 [42,]1 0 0 0 0 [43,]0 0 0 1 0 [44,]0 0 1 0 0 [45,]1 0 0 0 0 [46,]1 0 0 0 0 [47,]0 0 0 1 0 [48,]1 0 0 0 0 [49,]0 0 0 1 0 [50,]0 0 1 0 0 [51,]0 0 0 1 0 [52,]1 0 0 0 0 [53,]0 0 0 1 0 [54,]1 0 0 0 0 [55,]0 1 0 0 0 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Weibull survival modeling with covariate
I was wondering if someone familiar with survival analysis can help me with the following. I would like to fit a Weibull curve, that may be dependent on a covariate, my dataframe labdata that has the fields cov, time, and censor. Do I do the following? wieb-survreg(Surv(labdata$time, labadata$censor)~labdata$cov, dist=weibull) This returns: weib Call: survreg(formula = Surv(labdata$time, labdata$censor) ~ labdata$cov, dist = weibull) Coefficients: (Intercept) labdata$cov 8.091955112 0.001552897 Scale= 0.7532474 Loglik(model)= -12633.6 Loglik(intercept only)= -12734.8 Chisq= 202.41 on 1 degrees of freedom, p= 0 n= 5496 I am not quite sure how to use the output. I see that it gives the Scale parameter. How do I find the Shape paramater as a function of the covariate? Thank you, Steven --- - Steven Shechter PhD Candidate in Industrial Engineering University of Pittsburgh www.pitt.edu/~sms13 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [ESS] Re: [R] Strange characters in 2.1.0?
On Wed, 08-Jun-2005 at 11:20AM +0200, Martin Maechler wrote: | PaCo == Patrick Connolly [EMAIL PROTECTED] | on Wed, 8 Jun 2005 11:31:44 +1200 writes: | | PaCo On Tue, 07-Jun-2005 at 04:10PM +0200, Martin Maechler wrote: | PaCo | Dan == Dan Bolser [EMAIL PROTECTED] | | PaCo | .. | PaCo | | PaCo | Dan I have gone back to 2.0.0 :) | PaCo | | PaCo | Don't do that! | PaCo | You've lost tons of nice new features and gained quite an amount | PaCo | of old bugs by downgrading .. | | PaCo I get the non-generic quotes to show on the screen, but they won't | PaCo print with enscript. I end up with a lot of wrapped lines and | PaCo nonsense where an unknown character should be. | | Why is this diverted from R- to ESS-help? I erroneously thought that I could see the quotes only within Emacs and not in a terminal window, so I thought it was half way there with Emacs but not started in a terminal. I'll check more carefully henceforth. [...] | If I were in New Zeeland and would not need accents or umlauts, | I'd probably stick with latin1 (and would make sure my X | server got proper non-utf8 fonts) for another year or so. Thanks for the clarification. -- Patrick Connolly HortResearch Mt Albert Auckland New Zealand Ph: +64-9 815 4200 x 7188 ~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~ I have the world`s largest collection of seashells. I keep it on all the beaches of the world ... Perhaps you`ve seen it. ---Steven Wright ~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Forecasting with macroeconomic structural equations models?
Hello, Is there a package or sample code that shows how to do ex ante forecasts with a macroeconomic structural equations model? I looked at the sem package, which lets you estimate e.g. Klein's model, but I'm not sure how to make simulations using the full set of equations, including the identities. Thank you, Ronaldo Carpio [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Weibull survival modeling with covariate
Hi, I'm also wondering which expression the survreg() uses for Weibull regression. Referring to help(survreg) and help(survreg.distributions), I guess survreg() fits the following model. survreg() uses a different parametrization, say F(x, Wshape, Wscale) = 1-exp(-Wscale*(x^Wshape))), and fits a parametric model with these formulas. Wshape = 1/Scale (calculated by survreg()) log(Wscale) = model with covariates Is it correct? Thanks a lot. Watalu I was wondering if someone familiar with survival analysis can help me with the following. I would like to fit a Weibull curve, that may be dependent on a covariate, my dataframe labdata that has the fields cov, time, and censor. Do I do the following? wieb-survreg(Surv(labdata$time, labadata$censor)~labdata$cov, dist=weibull) This returns: weib Call: survreg(formula = Surv(labdata$time, labdata$censor) ~ labdata$cov, dist = weibull) Coefficients: (Intercept) labdata$cov 8.091955112 0.001552897 Scale= 0.7532474 Loglik(model)= -12633.6 Loglik(intercept only)= -12734.8 Chisq= 202.41 on 1 degrees of freedom, p= 0 n= 5496 I am not quite sure how to use the output. I see that it gives the Scale parameter. How do I find the Shape paramater as a function of the covariate? Thank you, Steven --- - Steven Shechter PhD Candidate in Industrial Engineering University of Pittsburgh www.pitt.edu/~sms13 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] How to plot more than 3 sets in Venn Diagrams?
If the data you posted is prototypical of your datasets then note that: - it has two disconnected components of 1 and 4 sets - there are only 8 unique rows out of 32 - 5 of these 8 are the regions that contain only the non-intersecting portion of each of the 5 sets unique(mydata) AA CH EA IN MY 1 0 0 0 1 0 2 1 0 0 0 0 4 0 0 0 0 1 8 1 0 0 1 0 14 0 0 1 0 0 21 0 1 0 0 0 24 1 0 1 1 1 37 0 0 1 1 0 Thus these are not very interesting Venn diagrams as they are quite sparse. Perhaps you would be better off to represent them as bipartite graphs with a node for each row and a node for each column with the edges corresponding to the ones. The RGraphViz package (or graphViz, itself, locatable via google) could be used for that. On 6/8/05, Tan Hui Hui Jenny [EMAIL PROTECTED] wrote: I'm trying to plot Venn diagrams with more than 3 sets (5 actually) in order to describe graphically the genetic variation between populations. I tried the limma library but realised it can only plot 3 sets. Is there any solution? Of course I could plot the chart manually but it'll take too long (have other datasets). One of my dataset is given below. THanks for any advice. j AA CH EA IN MY [1,] 0 0 0 1 0 [2,] 1 0 0 0 0 [3,] 1 0 0 0 0 [4,] 0 0 0 0 1 [5,] 1 0 0 0 0 [6,] 1 0 0 0 0 [7,] 1 0 0 0 0 [8,] 1 0 0 1 0 [9,] 1 0 0 0 0 [10,]1 0 0 0 0 [11,]1 0 0 0 0 [12,]0 0 0 0 1 [13,]1 0 0 0 0 [14,]0 0 1 0 0 [15,]1 0 0 0 0 [16,]0 0 1 0 0 [17,]1 0 0 0 0 [18,]0 0 0 1 0 [19,]1 0 0 1 0 [20,]0 0 1 0 0 [21,]0 1 0 0 0 [22,]1 0 0 0 0 [23,]0 0 1 0 0 [24,]1 0 1 1 1 [25,]0 1 0 0 0 [26,]1 0 0 0 0 [27,]1 0 0 0 0 [28,]0 0 0 1 0 [29,]0 0 0 1 0 [30,]1 0 0 0 0 [31,]1 0 0 0 0 [32,]0 0 0 1 0 [33,]0 0 0 1 0 [34,]0 1 0 0 0 [35,]1 0 0 0 0 [36,]0 0 0 1 0 [37,]0 0 1 1 0 [38,]1 0 0 1 0 [39,]0 0 0 1 0 [40,]0 0 0 1 0 [41,]0 1 0 0 0 [42,]1 0 0 0 0 [43,]0 0 0 1 0 [44,]0 0 1 0 0 [45,]1 0 0 0 0 [46,]1 0 0 0 0 [47,]0 0 0 1 0 [48,]1 0 0 0 0 [49,]0 0 0 1 0 [50,]0 0 1 0 0 [51,]0 0 0 1 0 [52,]1 0 0 0 0 [53,]0 0 0 1 0 [54,]1 0 0 0 0 [55,]0 1 0 0 0 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] May I ask you a question about matrix population models?
Have you found the file Using Rmetasim? In windows you can access this file by using the help and selecting browse directory. There appears to be a reasonable amount of information here. It looks to me as if you need to work your way through these files until you understand what is going on. Once you have done that you might be able to be more specific with your question. Tom -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] Behalf Of luan_sheng Sent: Thursday, 9 June 2005 10:27 AM To: R-help Subject: [R] May I ask you a question about matrix population models? Dear R user: Now I am studying matrix population models and Rmetasim package,but I find it's very difficult to understand this model fully and can't find a good teacher.My question is that:Now I do some research about one marine shrimp: Penaeus chinensis, it's life cycle is one year. after it breed it's offsprings, it will died. This shrimp's life cycle is dicrete ,if i use project matrix or life cyle graph to express it,how can i do it? Best regards luan_sheng __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
