Re: [R] Odp: Error in matrix, not ordered vectors or numerical value, and SIAR.

2011-11-11 Thread Alex DJ
Thank-you Petr. 

I have consulted as many manuals and help pages, and search engines, and
trying various things in SIAR, but continuous errors prevail.. 

I have tried changing all matrices to numeric, but the first column , of the
two latter matrices, returns as NA. Any other suggestions? 

I had put in one of these messages more detail about data, see below. 


pond - c(A3, A3, A3, B3, B3, B3,C2, A7, A7)
d13C - c(-22.07, -20.39, -21.11...) 
d15N - c(16.06, 17.13, 16.7, ...)
FConsumerData - data.frame(pond, d13C, d15N)
FTEFData - data.frame(tefsources, tefMeand13C, tefSDd13C, tefMeand15N,
tefSDd15N)
sources - c(Copepoda, Algal Paste, Amphipoda,
Meand13C - c(-26.57, -20.41, -21.65,...
SDd13C - c(1.78, 4.01, 1.46,...
Meand15N - c(15.55, 9.58, 14.07...
SDd15N - c(1.77, 1.34, 1.66,...
FSourceData - data.frame(sources, Meand13C, SDd13C, Meand15N, SDd15N)
tefsources - c(Copepoda, Algal Paste, Amphipoda,...
tefMeand13C - c(0.4, 0.4, 0.4,...
tefSDd13C - c(1.3, 1.3, 1.3,
tefMeand15N - c(3.4, 3.4, 3.4, 
tefSDd15N - c(1.0, 1.0, 1.0, 

Hope that helps explain? 

?? 

--
View this message in context: 
http://r.789695.n4.nabble.com/Error-in-matrix-not-ordered-vectors-or-numerical-value-and-SIAR-tp4024578p4030650.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Errors in SIAR

2011-11-11 Thread Mark Difford
On Nov 11, 2011 at 1:06am Alex DJ wrote:

Alex,

I haven't followed this thread closely (now split), but this most recent
request for help from you is very difficult to make sense of. I think we
need access to your data set. Further, there appear to be problems with the
function you are using in package siar[1], based on trying to run demos in
the package. (Though this may not be the issue.)

Two things you can try are the following: First, the function you are using
(siarmcmcdirichletv4) calls for matrices, not data frames (which you are
using). Usually this throws an error, so maybe this is not the problem.

Secondly, if you have a factor then you turn it into an ordered factor by
doing

ordered(myFactor)

where myFactor is the name of the factor you want to convert to an ordered
factor. You can call the levels of your factor what you like. Here, reading
up on these functions, and trying the examples, would help you to help
yourself. Have you even looked at ?factor. The fact that you are new to R
really means that you should be reading everything about R that you can lay
your hands on. This will help you to ask sensible questions about things
that still puzzle you.

And if you want to convert a factor into a numeric then you have to do
something like

as.numeric(as.character(myFactor))

This is all documented in R. Have you read the manualAn Introduction to R,
which comes with R? There is also plenty of good documentation on how to use
R on the web and indeed on R`s homepage.

Regards, Mark.

[1] R is case-sensitive: the package is called siar, not SIAR. Please
respect the author's designation.

-
Mark Difford (Ph.D.)
Research Associate
Botany Department
Nelson Mandela Metropolitan University
Port Elizabeth, South Africa
--
View this message in context: 
http://r.789695.n4.nabble.com/Errors-in-SIAR-tp4029804p4030667.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] Have write.csv write csv without observation number

2011-11-11 Thread dlarmour
Hello all,

Have a little annoying issue with wreite.csv. It always writes out
observation number. Is there a way to avoid this. Thanks,

DL

--
View this message in context: 
http://r.789695.n4.nabble.com/Have-write-csv-write-csv-without-observation-number-tp4030672p4030672.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Generating the Ctrl-M character

2011-11-11 Thread Ashim Kapoor
The \r\n worked as the the diff can see no difference between my program's
file and the file generated from the windows computer.

For David,

I use Emacs whitespace mode to see the spaces and the \n's. The \r does not
show up in the whitespace mode. Maybe there is some way of turning it on.

but THANK you guys.

Ashim : )

On Fri, Nov 11, 2011 at 12:11 PM, Jeff Newmiller
jdnew...@dcn.davis.ca.uswrote:

 No need to get defensive... I couldn't tell whether you remembered the
 details. This is starting to go OT, though.

 As to handing \n as expected, in the context of generating text files
 for their native platforms they do act as expected.  However, in the
 context of people writing binary characters in a file in an attempt to make
 it suitable for some other platform, munging can happen, and people can be
 surprised by the results if they don't get the whole picture.

 Specifically, AFAIK \n is always stored as a single character (LF) in
 memory on UTF8 or ASCII systems. On Windows, when that character is written
 to a text file a \r is written just before it, creating the expected
 end-of-line marker for a text file on that OS.  On the old Macintosh
 operating system a  lone \r was the EOL marker, so \n was automatically
 translated to \r (but not on OSX, which uses UNIX conventions). From the
 perspective of the person working with native text files, this makes things
 just work. From the perspective of someone trying to obtain a particular
 but pattern in the output file this amounts to munging.

 The C programming standard library includes the concept of distinguishing
 between text files and binary files just so this native adaptation can
 occur or not, depending on the desired behavior.
 ---
 Jeff NewmillerThe .   .  Go Live...
 DCN:jdnew...@dcn.davis.ca.usBasics: ##.#.   ##.#.  Live
 Go...
  Live:   OO#.. Dead: OO#..  Playing
 Research Engineer (Solar/BatteriesO.O#.   #.O#.  with
 /Software/Embedded Controllers)   .OO#.   .OO#.  rocks...1k
 ---
 Sent from my phone. Please excuse my brevity.

 David Winsemius dwinsem...@comcast.net wrote:

 
 On Nov 10, 2011, at 11:55 PM, Jeff Newmiller wrote:
 
  Wow, Deadpan David.
 
  How about using the escape sequence \r?
 
 So shoot me. It wasn't documented in a manner that I recognized in the
 
 places I looked.
 
 ?character
 ?Syntax
 ?Constants
 
 And I did look at ?Quotes where \r is listed but did not know that
 it was == cntrl-M (if in fact it is.) I assumed (probably incorrectly)
 
 that \r was cntrl-R. There was proabably a time in the past when I
 could have told you a lot of the decimal and maybe even hexadecimal
 equivalents for cr, lf, beep, but those days are behind me.
 
 
  Keep in mind that Ctrl-M is used as the end-of-line character on
  some operating systems, so accomplishing this may not be portable,
  and you didn't specify your operating system. On the three main
  platforms (*nix, Windows/DOS, and Mac), \r should work, but \n
  may get munged.
 
 What? I very much doubt that any of those systems will not handle \n
 
 as expected. \r on the other hand I'm not so sure of.
 
 --
 David.
 
 
 

 ---
  Jeff NewmillerThe .   .  Go
  Live...
  DCN:jdnew...@dcn.davis.ca.usBasics: ##.#.   ##.#.
  Live Go...
   Live:   OO#.. Dead: OO#..
  Playing
  Research Engineer (Solar/BatteriesO.O#.   #.O#.  with
  /Software/Embedded Controllers)   .OO#.   .OO#.
  rocks...1k
 

 ---
  Sent from my phone. Please excuse my brevity.
 
  David Winsemius dwinsem...@comcast.net wrote:
 
 
  On Nov 10, 2011, at 9:35 PM, Ashim Kapoor wrote:
 
  Dear R-helpers,
 
  I want to append a Ctrl-M character to a string and then save it to
  a text
  file.
 
  mystring-This is a test.
 
  # How do I add a Ctrl-M to it in the end ??
 
  cat(mystring,file=testfile)
 
 
  cntrl_m - intToUtf8(13)
 
  cat(cntrl_m,file=testfile)
 
  The resulting file seems to have a blank line in my editor.
 
  --
 
  David Winsemius, MD
  West Hartford, CT
 
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide
  http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
 
 
 David Winsemius, MD
 West Hartford, CT



[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting 

Re: [R] Have write.csv write csv without observation number

2011-11-11 Thread Timothy Bates
If you find yourself asking is there a way..?, put that question mark in 
front of the function you are wondering about and push return

?write.csv

In this case, you will see that the write function has a parameter to say 
whether to write out row names ornot, which defaults to TRUE...

 row.names = TRUE


On Nov 11, 2011, at 7:46 AM, dlarmour wrote:

 Hello all,
 
 Have a little annoying issue with wreite.csv. It always writes out
 observation number. Is there a way to avoid this. Thanks,
 
 DL
 
 --
 View this message in context: 
 http://r.789695.n4.nabble.com/Have-write-csv-write-csv-without-observation-number-tp4030672p4030672.html
 Sent from the R help mailing list archive at Nabble.com.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Error with Source()

2011-11-11 Thread Duncan Murdoch

On 11-11-10 9:25 AM, ftonini wrote:

Hi everybody,

I started to receive a weird message in R that I have never seen
before...also I haven't found anything on google or on this forum about it.
Whenever I use the command source(...) to point to one of my scripts, I get
the following message:

Error in source(myfunctions.R) : myfunctions.R:884:9: unexpected symbol
883:
884: cond
  ^

I am using the same commands as I did in the past and it was working...I
started to receive this error (not sure if it has to do with it or not)
after trying to create a batch file to run one of my .R scripts with
double-click. That batch file worked...but as soon as I use the source()
command it does not work any more.

Any help is appreciated!


The message indicates that on line 884 of that script, the parser is 
finding a syntax error.  A common way to generate that error is to fail 
to close parentheses, e.g.


( a
cond

gives a similar one.

Duncan Murdoch

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Gibbs sampler

2011-11-11 Thread Duncan Murdoch

On 11-11-10 1:28 PM, Gyanendra Pokharel wrote:

I have the following code,
gibbs-function(m,theta = 0.25, lambda =0.55, n =1){
 alpha- 1.5
 beta- 1.5
 gamma- 1.5
 x- array(0,c(m+1, 3))
 x[1,1]- theta
 x[1,2]- lambda
 x[1,3]- n
 for(t in 2:(m+1)){
 x[t,1]- rbinom(1, x[t-1,3], x[t-1,1])
 x[t,2]-rbeta(1, x[t-1,1] + alpha, x[t-1,3] - x[t-1,1] + beta)


Are you certain that x[t-1,3] - x[t-1,1] + beta will always be positive? 
 If it is negative, rbeta will return NaN.


Duncan Murdoch


 x[t,3]- rpois(1,(1 - x[t-1,1])*gamma)
 }
 x
}
gibbs(100)

it returns only 1 or 2 values of theta, some times NaN, this may be if any
theta is greater than 1, which is used as the probability for the next
rbinom(), so returns NaN. Can some one suggest to solve this problem?

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] beanplot without log scale?

2011-11-11 Thread Dennis Murphy
Hi:

Try the log =  argument to beanplot(). Here's an example:

library('beanplot')
v - exp(runif(1000, 0, 10))
beanplot(v)
beanplot(v, log = ')

HTH,
Dennis

On Thu, Nov 10, 2011 at 1:27 PM, Twila Moon twi...@uw.edu wrote:
 Is it possible to force beanplot not to use a log scale? I want to be able to 
 create multiple beanplots of different data on the same specified y-axis.

 When I try to force a ylim, I get the following...
 beanplot(reg5vel,ylim=(c(0,12000)))
 log=y selected
 Warning message:
 In plot.window(xlim = c(0.5, 7.5), ylim = c(0, 12000), log = y) :
  nonfinite axis limits [GScale(-inf,4.07918,2, .); log=1]

 Thanks,
 Twila

        [[alternative HTML version deleted]]

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] performance of adaptIntegrate vs. integrate

2011-11-11 Thread Hans W Borchers
baptiste auguie baptiste.auguie at googlemail.com writes:

 
 Dear list,
 
 [cross-posting from Stack Overflow where this question has remained
 unanswered for two weeks]
 
 I'd like to perform a numerical integration in one dimension,
 
 I = int_a^b f(x) dx
 
 where the integrand f: x in IR - f(x) in IR^p is vector-valued.
 integrate() only allows scalar integrands, thus I would need to call
 it many (p=200 typically) times, which sounds suboptimal. The cubature
 package seems well suited, as illustrated below,
 
 library(cubature)
 Nmax - 1e3
 tolerance - 1e-4
 integrand - function(x, a=0.01) c(exp(-x^2/a^2), cos(x))
 adaptIntegrate(integrand, -1, 1, tolerance, 2, max=Nmax)$integral
 [1] 0.01772454 1.68294197
 
 However, adaptIntegrate appears to perform quite poorly when compared
 to integrate. Consider the following example (one-dimensional
 integrand),
 
 library(cubature)
 integrand - function(x, a=0.01) exp(-x^2/a^2)*cos(x)
 Nmax - 1e3
 tolerance - 1e-4
 
 # using cubature's adaptIntegrate
 time1 - system.time(replicate(1e3, {
   a - adaptIntegrate(integrand, -1, 1, tolerance, 1, max=Nmax)
 }) )
 
 # using integrate
 time2 - system.time(replicate(1e3, {
   b - integrate(integrand, -1, 1, rel.tol=tolerance, subdivisions=Nmax)
 }) )
 
 time1
 user  system elapsed
   2.398   0.004   2.403
 time2
 user  system elapsed
   0.204   0.004   0.208
 
 a$integral
  [1] 0.0177241
 b$value
  [1] 0.0177241
 
 a$functionEvaluations
  [1] 345
 b$subdivisions
  [1] 10
 
 Somehow, adaptIntegrate was using many more function evaluations for a
 similar precision. Both methods apparently use Gauss-Kronrod
 quadrature, though ?integrate adds a Wynn's Epsilon algorithm. Could
 that explain the large timing difference?


Cubature is astonishingly slow here though it was robust and accurate in
most cases I used it. You may write to the maintainer for more information.

Even a pure R implementation of adaptive Gauss-Kronrod as in pracma::quadgk
is faster than cubature:

library(pracma)
time3 - system.time(replicate(1e3, 
  { d - quadgk(integrand, -1, 1) }# 0.0177240954011303
) )
#  user  system  elapsed 
# 0.899   0.0020.897

If your functions are somehow similar and you are willing to dispense with
the adaptive part, then compute Gaussian quadrature nodes xi and weights wi
for the fixed interval and compute sum(wi * fj(xi)) for each function fj.
This avoids recomputing nodes and weights for every call or function.

Package 'gaussquad' provides such nodes and weights. Or copy the relevant
calculations from quadgk (the usual (7,15)-rule).

Regards, Hans Werner


 I'm open to suggestions of alternative ways of dealing with
 vector-valued integrands.
 
 Thanks.
 
 baptiste
 


__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Building a statically-linked extension?

2011-11-11 Thread Prof Brian Ripley

On Thu, 10 Nov 2011, Tyler Pirtle wrote:


Hi there,

I'm writing an R extension that has a C component that relies on two third
party libraries that I'm bundling
with the extension.

I'd like to statically link my resulting extension so that I can rely on
the bundled versions of the libraries I'm
distributing, so there's two questions -

1, does R allow statically linked C extensions to be used at runtime?


Yes


2, are there any standard ways of having R build my extension statically?


Yes.

If you want to know more, follow the posting guide: this is not an 
R-help topic, the missing 'at a minimum' information is vital 




Thanks,


Tyler

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.



--
Brian D. Ripley,  rip...@stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] match first consecutive list of capitalized words in string

2011-11-11 Thread Richter-Dumke, Jonas
Thank you very much for your help.

I'm using the second suggestion in my program and it works very well.

Jonas


-Original Message-
From: Gabor Grothendieck [mailto:ggrothendi...@gmail.com]
Sent: Thu 11/10/2011 5:58 AM
To: Richter-Dumke, Jonas
Cc: r-help@r-project.org
Subject: Re: [R] match first consecutive list of capitalized words in string
 
On Tue, Nov 8, 2011 at 7:48 AM, Richter-Dumke, Jonas
rich...@demogr.mpg.de wrote:
 Dear R-Helpers,

 this is my first post ever to a mailing list, so please feel free to point 
 out any missunderstandings on my side regarding the conventions of this 
 mailing list.

 My problem:

 Assuming the following character vector is given:

 names - c(filia Maria, vidua Joh Dirck Kleve (oo 02.02.1732), Bernardus 
 Engelb Franciscus Linde j.u.Doktor referendarius sereniss Judex et gograven 
 Rheinensis)

 Is there a regular expression matching the first consecutive list of 
 capitalized words in a single characterstring (Maria, Joh Dirck Kleve, 
 Bernardus Engelb Franciscus Linde)?
 This expression would very reliably seperate the person names from the 
 additional information in my historic church register transcription.


Try this. It matches a word boundary followed by zero or more of the
parenthesized expression.  That expression is an upper case letter
followed by zero or more lower case letters followed by one or more
spaces.  Finally we match the last word which consists of an upper
case letter followed by zero or more lower case letters and a word
boundary.  Note that it assumes R 2.14.0 or later:

 re - \\b([[:upper:]][[:lower:]]* +)*[[:upper:]][[:lower:]]*\\b
 regmatches(names, regexpr(re, names))
[1] Maria Joh Dirck Kleve
[3] Bernardus Engelb Franciscus Linde

-- 
Statistics  Software Consulting
GKX Group, GKX Associates Inc.
tel: 1-877-GKX-GROUP
email: ggrothendieck at gmail.com


--
This mail has been sent through the MPI for Demographic ...{{dropped:10}}

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] Win upgrade pb (virus)

2011-11-11 Thread Mario Valle
I just upgraded my Win7 32bits installation to 2.14.0 after deinstalling 
2.12.x
First thing I moved the win-library from 2.12 to 2.14 and executed a 
update.packages(ask='graphics',checkBuilt=TRUE)

(Swiss mirror).
This aborts with the console message:

Error in if (any(diff)) { : missing value where TRUE/FALSE needed

And the antivirus (AVG) pops up complaining that colorspace.dll contains 
the virus Win32/Heur
I cannot ignore the thread because update.packages has already crashed 
when I can reach the antivirus dialogbox to push the ignore button.

To continue I had to remove the colorspace and ggplot2 packages.

Is it a known problem? What else can I do (except removing the above 
packages or changing antivirus)?


Thanks!
mario

--
Ing. Mario Valle
Data Analysis and Visualization Group| http://www.cscs.ch/~mvalle
Swiss National Supercomputing Centre (CSCS)  | Tel:  +41 (91) 610.82.60
v. Cantonale Galleria 2, 6928 Manno, Switzerland | Fax:  +41 (91) 610.82.82

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Win upgrade pb (virus)

2011-11-11 Thread Prof Brian Ripley

On Fri, 11 Nov 2011, Mario Valle wrote:

I just upgraded my Win7 32bits installation to 2.14.0 after deinstalling 
2.12.x
First thing I moved the win-library from 2.12 to 2.14 and executed a 
update.packages(ask='graphics',checkBuilt=TRUE)

(Swiss mirror).
This aborts with the console message:

Error in if (any(diff)) { : missing value where TRUE/FALSE needed

And the antivirus (AVG) pops up complaining that colorspace.dll contains the 
virus Win32/Heur
I cannot ignore the thread because update.packages has already crashed when I 
can reach the antivirus dialogbox to push the ignore button.

To continue I had to remove the colorspace and ggplot2 packages.

Is it a known problem? What else can I do (except removing the above packages 
or changing antivirus)?


Yes (a known problem in AVG), and the standard advice is to change 
antivirus (or set exceptions, which I gather AVG does not allow). 
And please report the false positive to AVG: the more they hear about 
the more attention they will give it.


The binary distribution has been checked by e.g. Sophos (which both 
winbuilder and ox.ac.uk use).




Thanks!
   mario

--
Ing. Mario Valle
Data Analysis and Visualization Group| http://www.cscs.ch/~mvalle
Swiss National Supercomputing Centre (CSCS)  | Tel:  +41 (91) 610.82.60
v. Cantonale Galleria 2, 6928 Manno, Switzerland | Fax:  +41 (91) 610.82.82



--
Brian D. Ripley,  rip...@stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] optim seems to be finding a local minimum

2011-11-11 Thread Hans W Borchers
Ben Bolker bbolker at gmail.com writes:

 
   Simulated annealing and other stochastic global optimization 
 methods are also possible solutions, although they may or may not
 work better than the many-starting-points solution -- it depends
 on the problem, and pretty much everything has to be tuned.  Tabu
 search http://en.wikipedia.org/wiki/Tabu_search is another possibility,
 although I don't know much about it ...
 

It is known that the Excel Solver has much improved during recent years.
Still there are slightly better points such as

myfunc(c(0.889764228112319, 94701144.5712312))   # 334.18844

restricting the domain to [0, 1] x [0, 10^9] for an evolutionary approach,
for instance DEoptim::DEoptim().

Finding a global optimum in 2 dimensions is not so difficult. Here the scale
of the second variable could pose a problem as small local minima might be
overlooked easily.

Hans Werner

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] R package for segmentation with both continuous and categorical input variables XXXX

2011-11-11 Thread Ingmar Visser
Th CRAN taskview on clustering has many clustering and mixture packages
which
may be useful for your precise situation:

http://cran.r-project.org/web/views/Cluster.html

hth, Ingmar

On Thu, Nov 10, 2011 at 7:41 PM, Dan Abner dan.abne...@gmail.com wrote:

 Hello everyone,

 Can anyone suggest a decently documented (with good examples in the
 documentation) R package/function that performs segmentation (cluster,
 mixture modeling) of a population using both continuous and categorical
 input variables?

 Thank you,

 Dan

[[alternative HTML version deleted]]

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Odp: Error in matrix, not ordered vectors or numerical value, and SIAR.

2011-11-11 Thread Petr PIKAL
 
 Re: [R] Odp:  Error in matrix, not ordered vectors or numerical value, 
and SIAR.
 
 Thank-you Petr. 
 
 I have consulted as many manuals and help pages, and search engines, and
 trying various things in SIAR, but continuous errors prevail.. 
 
 I have tried changing all matrices to numeric, but the first column , of 
the
 two latter matrices, returns as NA. Any other suggestions? 
 
 I had put in one of these messages more detail about data, see below. 
 
 
 pond - c(A3, A3, A3, B3, B3, B3,C2, A7, A7)
 d13C - c(-22.07, -20.39, -21.11...) 
 d15N - c(16.06, 17.13, 16.7, ...)
 FConsumerData - data.frame(pond, d13C, d15N)
 FTEFData - data.frame(tefsources, tefMeand13C, tefSDd13C, tefMeand15N,
 tefSDd15N)
 sources - c(Copepoda, Algal Paste, Amphipoda,
 Meand13C - c(-26.57, -20.41, -21.65,...
 SDd13C - c(1.78, 4.01, 1.46,...
 Meand15N - c(15.55, 9.58, 14.07...
 SDd15N - c(1.77, 1.34, 1.66,...
 FSourceData - data.frame(sources, Meand13C, SDd13C, Meand15N, SDd15N)
 tefsources - c(Copepoda, Algal Paste, Amphipoda,...
 tefMeand13C - c(0.4, 0.4, 0.4,...
 tefSDd13C - c(1.3, 1.3, 1.3,
 tefMeand15N - c(3.4, 3.4, 3.4, 
 tefSDd15N - c(1.0, 1.0, 1.0, 
 
 Hope that helps explain? 

Not much. 

FConsumerData is data frame with nonnumeric first column, the same is with 
FTEFData and FSourceData. If you change it to matrix you will get 
character matrix. Trying to change character values to numeric obviously 
leads to NA. How the poor R shall know what number it shall make from A3, 
B3, C2 or A7?

I do not use siar package but in the help page it is clearly stated that 
input shall be ***numeric*** matrices. Not data frames, not character 
values not factors. If you do not know distinction among those objects you 
shall look into R intro document (at least first 6 chapters) which should 
came along with your installation.

BTW, in help page there is

out - 
siarmcmcdirichletv4(geese1demo,sourcesdemo,correctionsdemo,concdepdemo)

It does not produce any error does it?

If it does not produce any error it is sign that something is wrong with 
your data. Look at structure of your objects you are trying to evaluate 
and compare it with those demo objects.

str(object)

is the first option when something does not seem to be OK.

Regards
Petr




 
 ?? 
 
 --
 View this message in context: http://r.789695.n4.nabble.com/Error-in-
 
matrix-not-ordered-vectors-or-numerical-value-and-SIAR-tp4024578p4030650.html
 Sent from the R help mailing list archive at Nabble.com.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide 
http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Building package problem

2011-11-11 Thread Eduardo Mendes
Hello

Many thanks for the replies.

I am note sure whether you've got what you meant (Prof. Ripley) but here is
the output of Sys.getlocale()

 Sys.getlocale()[1] LC_COLLATE=English_United 
 States.1252;LC_CTYPE=English_United States.1252;LC_MONETARY=English_United 
 States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252


Is that what you meant?

Many thanks

Ed


On Tue, Nov 8, 2011 at 1:03 PM, Prof Brian Ripley rip...@stats.ox.ac.ukwrote:

 R CMD check is *not* 'building a package'.  Nor is making a Windows binary
 package.  'Building a package' is creating a source tarball from a source
 directory.


 On Tue, 8 Nov 2011, Joshua Wiley wrote:

  Hi Ed,

 If the only error is in examples then this should work:

 R CMD check --no-examples foopkg

 should not have anything to do with vignettes (although those may also
 not run, who knows).  As far as building a binary, look at:

 R CMD INSTALL --help

 which leads you to

 R CMD INSTALL --build foopkg


 And as for the hydroGOF issue, my guess is that the problem is your locale
 or timezone.  But despite the posting guide, you failed to tell us. AFAIK
 CRAN only checks packages in English locales.

 (One thing we know is that in Columbia there was no midnight on one of the
 dates in that file.  So hydroGOF really ought to be specifying a timezone
 when reading character datetimes.)



 HTH,

 Josh



 On Tue, Nov 8, 2011 at 4:35 AM, Eduardo M. A. M. Mendes
 emammen...@gmail.com wrote:

 Dear R-users

 I am trying to recompile a CRAN package on Windows 32.   Rtools for 2.14
 (that is the version I am running) and miktex were sucessfully installed on
 my machine.

 Problems:

 a) hydroGOF is a CRAN package, but R CMD check does not work on it.

 C:\Users\eduardo\Documents\R_**tests2R CMD check hydroGOF
 * using log directory 'C:/Users/eduardo/Documents/R_**
 tests2/hydroGOF.Rcheck'
 * using R version 2.14.0 (2011-10-31)
 * using platform: i386-pc-mingw32 (32-bit)
 * using session charset: ISO8859-1
 * checking for file 'hydroGOF/DESCRIPTION' ... OK
 * checking extension type ... Package
 * this is package 'hydroGOF' version '0.3-2'
 * checking package namespace information ... OK
 * checking package dependencies ... OK
 * checking if this is a source package ... OK
 * checking if there is a namespace ... OK
 * checking for executable files ... OK
 * checking whether package 'hydroGOF' can be installed ... OK
 * checking installed package size ... OK
 * checking package directory ... OK
 * checking for portable file names ... OK
 * checking DESCRIPTION meta-information ... OK
 * checking top-level files ... OK
 * checking index information ... OK
 * checking package subdirectories ... OK
 * checking R files for non-ASCII characters ... OK
 * checking R files for syntax errors ... OK
 * checking whether the package can be loaded ... OK
 * checking whether the package can be loaded with stated dependencies
 ... OK
 * checking whether the package can be unloaded cleanly ... OK
 * checking whether the namespace can be loaded with stated dependencies
 ... OK
 * checking whether the namespace can be unloaded cleanly ... OK
 * checking for unstated dependencies in R code ... OK
 * checking S3 generic/method consistency ... OK
 * checking replacement functions ... OK
 * checking foreign function calls ... OK
 * checking R code for possible problems ... OK
 * checking Rd files ... OK
 * checking Rd metadata ... OK
 * checking Rd cross-references ... OK
 * checking for missing documentation entries ... OK
 * checking for code/documentation mismatches ... OK
 * checking Rd \usage sections ... OK
 * checking Rd contents ... OK
 * checking for unstated dependencies in examples ... OK
 * checking contents of 'data' directory ... OK
 * checking data for non-ASCII characters ... OK
 * checking data for ASCII and uncompressed saves ... OK
 * checking examples ... ERROR
 Running examples in 'hydroGOF-Ex.R' failed
 The error most likely occurred in:

  ### Name: plot2
 ### Title: Plotting 2 Time Series
 ### Aliases: plot2
 ### Keywords: dplot

 ### ** Examples

 sim - 2:11
 obs - 1:10
 ## Not run:
 ##D plot2(sim, obs)
 ## End(Not run)

 ##
 # Loading daily streamflows of the Ega River (Spain), from 1961 to 1970
 require(zoo)

 Loading required package: zoo

 Attaching package: 'zoo'

 The following object(s) are masked from 'package:base':

as.Date, as.Date.numeric

  data(EgaEnEstellaQts)
 obs - EgaEnEstellaQts

 # Generating a simulated daily time series, initially equal to the
 observed se

 ries

 sim - obs

 # Randomly changing the first 2000 elements of 'sim', by using a normal
 distri

 bution

 # with mean 10 and standard deviation equal to 1 (default of 'rnorm').
 sim[1:2000] - obs[1:2000] + rnorm(2000, mean=10)

 # Plotting 'sim' and 'obs' in 2 separate panels
 plot2(x=obs, y=sim)

 # Plotting 'sim' and 'obs' in the same window
 plot2(x=obs, y=sim, plot.type=single)

 Error in as.POSIXlt.character(x, tz, ...) :
  character string is not in a 

Re: [R] A question on Programming

2011-11-11 Thread Oliver
Christofer Bogaso bogaso.christofer at gmail.com writes:
[...]
 Here my question is, is there any speed reduction if I put them within a
 function (I think there may be some speed reduction at least within
 for-loop, because that loop needs to call that function many times),
 relative to if I used that group of codes as-it-is in many places?
[...]


You did not asked for it, you may know it, or may not know it:
if you use apply functions and other vector oriented functions,
this can bring you a huge speedup, compared to a for-loop.


Ciao,
   Oliver

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Building package problem

2011-11-11 Thread Prof Brian Ripley

The posting guide asks for the output of sessionInfo() 
And what does Sys.timezone() say (it isn't always helpful).


On Fri, 11 Nov 2011, Eduardo Mendes wrote:


Hello
Many thanks for the replies.

I am note sure whether you've got what you meant (Prof. Ripley) but here is the
output of Sys.getlocale()

 Sys.getlocale()
[1] LC_COLLATE=English_United States.1252;LC_CTYPE=English_United 
States.1252;LC
_MONETARY=English_United States.1252;LC_NUMERIC=C;LC_TIME=English_United 
States.1
252

Is that what you meant?

Many thanks

Ed


On Tue, Nov 8, 2011 at 1:03 PM, Prof Brian Ripley rip...@stats.ox.ac.uk wrote:
  R CMD check is *not* 'building a package'.  Nor is making a Windows
  binary package.  'Building a package' is creating a source tarball
  from a source directory.

  On Tue, 8 Nov 2011, Joshua Wiley wrote:

Hi Ed,

If the only error is in examples then this should work:

R CMD check --no-examples foopkg

should not have anything to do with vignettes (although
those may also
not run, who knows).  As far as building a binary, look
at:

R CMD INSTALL --help

which leads you to

R CMD INSTALL --build foopkg


And as for the hydroGOF issue, my guess is that the problem is your locale
or timezone.  But despite the posting guide, you failed to tell us. AFAIK
CRAN only checks packages in English locales.

(One thing we know is that in Columbia there was no midnight on one of the
dates in that file.  So hydroGOF really ought to be specifying a timezone
when reading character datetimes.)



  HTH,

  Josh



  On Tue, Nov 8, 2011 at 4:35 AM, Eduardo M. A. M. Mendes
  emammen...@gmail.com wrote:
Dear R-users

I am trying to recompile a CRAN package on Windows
32.   Rtools for 2.14 (that is the version I am
running) and miktex were sucessfully installed on
my machine.

Problems:

a) hydroGOF is a CRAN package, but R CMD check does
not work on it.

C:\Users\eduardo\Documents\R_tests2R CMD check
hydroGOF
* using log directory
'C:/Users/eduardo/Documents/R_tests2/hydroGOF.Rcheck'
* using R version 2.14.0 (2011-10-31)
* using platform: i386-pc-mingw32 (32-bit)
* using session charset: ISO8859-1
* checking for file 'hydroGOF/DESCRIPTION' ... OK
* checking extension type ... Package
* this is package 'hydroGOF' version '0.3-2'
* checking package namespace information ... OK
* checking package dependencies ... OK
* checking if this is a source package ... OK
* checking if there is a namespace ... OK
* checking for executable files ... OK
* checking whether package 'hydroGOF' can be
installed ... OK
* checking installed package size ... OK
* checking package directory ... OK
* checking for portable file names ... OK
* checking DESCRIPTION meta-information ... OK
* checking top-level files ... OK
* checking index information ... OK
* checking package subdirectories ... OK
* checking R files for non-ASCII characters ... OK
* checking R files for syntax errors ... OK
* checking whether the package can be loaded ... OK
* checking whether the package can be loaded with
stated dependencies ... OK
* checking whether the package can be unloaded
cleanly ... OK
* checking whether the namespace can be loaded with
stated dependencies ... OK
* checking whether the namespace can be unloaded
cleanly ... OK
* checking for unstated dependencies in R code ...
OK
* checking S3 generic/method consistency ... OK
* checking replacement functions ... OK
* checking foreign function calls ... OK
* checking R code for possible problems ... OK
* checking Rd files ... OK
* checking Rd metadata ... OK
* checking Rd cross-references ... OK
* checking for missing documentation entries ... OK
* checking for code/documentation mismatches ... OK
* checking Rd \usage sections ... OK
* checking Rd contents ... OK
* checking for unstated dependencies in examples
... OK
* checking contents of 'data' directory ... OK
* checking data for non-ASCII characters ... OK
* checking data for ASCII and uncompressed saves
... OK
* checking examples ... ERROR
Running examples in 'hydroGOF-Ex.R' failed
The error most likely occurred in:

 

Re: [R] optim seems to be finding a local minimum

2011-11-11 Thread Ben Bolker
Hans W Borchers hwborchers at googlemail.com writes:

 
 Ben Bolker bbolker at gmail.com writes:
 
  
Simulated annealing and other stochastic global optimization 
  methods are also possible solutions, although they may or may not
  work better than the many-starting-points solution -- it depends
  on the problem, and pretty much everything has to be tuned.  Tabu
  search http://en.wikipedia.org/wiki/Tabu_search is another possibility,
  although I don't know much about it ...
  
 
 It is known that the Excel Solver has much improved during recent years.
 Still there are slightly better points such as
 
 myfunc(c(0.889764228112319, 94701144.5712312))   # 334.18844
 
 restricting the domain to [0, 1] x [0, 10^9] for an evolutionary approach,
 for instance DEoptim::DEoptim().
 
 Finding a global optimum in 2 dimensions is not so difficult. Here the scale
 of the second variable could pose a problem as small local minima might be
 overlooked easily.
 

  Have taken a (quick) second look at the problem, I agree that scaling
and centering are more likely to be useful solutions than stochastic
global optimization stuff.  Even using 'parscale' (see the optim
help page) may clear up the problem.

  Ben Bolker

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Title for a group of plots?

2011-11-11 Thread Uwe Ligges


Or even simpler (unbelievable!), see ?title:


par(mfcol=c(2,2), oma=c(0,0,1,0))
replicate(4, plot(1))
title(Hello World!, outer=TRUE)

Best,
Uwe Ligges



On 10.11.2011 17:13, Sarah Goslee wrote:

Try ?mtext for information on how to plot into the outer margin of a group of
plots.

Sarah

On Thu, Nov 10, 2011 at 11:07 AM, Kevin Burtonrkevinbur...@charter.net  wrote:

I can get multiple plots on a page like:



op- par(mfcol = c(3, 1))



What I was wondering is if there is a way to have a title for the whole
page? I can specify the title for each individual plot like:



plot(xxx, main=.)



But I would like a 'title' for the group of plots. Is this possible?



Thank you.



Kevin







__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] optim seems to be finding a local minimum

2011-11-11 Thread John C Nash
I won't requote all the other msgs, but the latest (and possibly a bit glitchy) 
version of
optimx on R-forge

1) finds that some methods wander into domains where the user function fails 
try() (new
optimx runs try() around all function calls). This includes L-BFGS-B

2) reports that the scaling is such that you really might not expect to get a 
good solution

then

3) Actually gets a better result than the

 xlf-myfunc(c(0.888452533990788,94812732.0897449))
 xlf
[1] 334.607


with Kelley's variant of Nelder Mead (from dfoptim package), with

 myoptx
  methodpar   fvalues fns  grs itns conv  KKT1
4 LBFGSB NA, NA 8.988466e+307  NA NULL NULL NA
2 Rvmmin   0.1, 200186870.6  25593.83  201 NULL0 FALSE
3 bobyqa 6.987875e-01, 2.001869e+08  1933.229  44   NA NULL0 FALSE
1   nmkb 8.897590e-01, 9.470163e+07  334.1901 204   NA NULL0 FALSE
   KKT2 xtimes  meths
4NA   0.01 LBFGSB
2 FALSE   0.11 Rvmmin
3 FALSE   0.24 bobyqa
1 FALSE   1.08   nmkb

But do note the terrible scaling. Hardly surprising that this function does not 
work. I'll
have to delve deeper to see what the scaling setup should be because of the 
nature of the
function setup involving some of the data. (optimx includes parscale on all 
methods).

However, original poster DID include code, so it was easy to do a quick check. 
Good for him.

JN

 ## Comparing this solution to Excel Solver solution:
 myfunc(c(0.888452533990788,94812732.0897449))
 
 -- Dimitri Liakhovitski marketfusionanalytics.com

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Sum of the deviance explained by each term in a gam model does not equal to the deviance explained by the full model.

2011-11-11 Thread Simon Wood
It's the same problem as with any regression model where the predictors 
are not strictly orthogonal. You can make the explained deviances per 
term sum to the whole model explained deviance by dropping terms 
sequentially, but then the explained deviance per term depends on the 
order in which you drop terms. The alternative in the code you pasted in 
below is at least independent of the ordering.


I suppose you could make things add up in the way you want by 
symmetrising the computation like this...


## reduce model one way
dev.1- (deviance(b1)-deviance(b))/deviance(b0) ## prop expl by s(x2)
dev.2- (deviance(b0)-deviance(b1))/deviance(b0) ## prop expl by  s(x1)

## reduce it the other way...
dev.2[2]- (deviance(b2)-deviance(b))/deviance(b0) ## prop expl by s(x2)
dev.1[2]- (deviance(b0)-deviance(b2))/deviance(b0) # prop expl by s(x1)

## average the alternatives
dev.1 - mean(dev.1)
dev.2 - mean(dev.2)

## so now explained deviance adds up...
dev.1+dev.2
summary(b)$dev.expl

... but I'm not convinced that this really adds anything useful, and it 
gets very tedious as the number of model terms increases...


--
Simon Wood, Mathematical Science, University of Bath BA2 7AY UK
+44 (0)1225 386603   http://people.bath.ac.uk/sw283


On 10/11/11 14:02, Huidong TIAN wrote:

Dear R users,
 I read your methods of extracting the variance explained by each
predictor in different places. My question is: using the method you
suggested, the sum of the deviance explained by all terms is not equal to
the deviance explained by the full model. Could you tell me what caused
such problem?


  set.seed(0)
  n-400
  x1- runif(n, 0, 1)
  ## to see problem with not fixing smoothing parameters
  ## remove the `##' from the next line, and the `sp'
  ## arguments from the `gam' calls generating b1 and b2.
  x2- runif(n, 0, 1) ## *.1 + x1
  f1- function(x) exp(2 * x)
  f2- function(x) 0.2*x^11*(10*(1-x))^6+10*(10*x)^3*(1-x)^10
  f- f1(x1) + f2(x2)
  e- rnorm(n, 0, 2)
  y- f + e
  ## fit full and reduced models...
  b- gam(y~s(x1)+s(x2))
  b1- gam(y~s(x1),sp=b$sp[1])
  b2- gam(y~s(x2),sp=b$sp[2])
  b0- gam(y~1)
  ## calculate proportions deviance explained...
  dev.1- (deviance(b1)-deviance(b))/deviance(b0) ## prop explained by

s(x2)

  dev.2- (deviance(b2)-deviance(b))/deviance(b0) ## prop explained by

s(x1)


  dev.1 + dev.2

[1] 0.6974949

  summary(b)$dev.expl

[1] 0.7298136

I checked the two models (b1  b2), found the model coefficients are
different with model b, so I feel it could be the problem.

wish to hear your comments.

Huidong Tian

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.



__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] proportional area venn diagrams?

2011-11-11 Thread Michael Friendly


In http://finzi.psych.upenn.edu/R/Rhelp02a/archive/14637.html

some rudimentary R functions were given  for drawing
proportional area venn diagrams with area of each intersection ~ the 
count in a 2 x 2 x 2 table.


I'm interested in this, for another application: showing the 
correlations among Y, X1, X2 using

area ~ r^2 of each pair (sometimes called a Ballantine diagram).

Before I attempt to hack the code given in that post, are there any 
packages/functions that do venn diagrams
with proportional areas?  Rseek turned up quite a view venn-like 
functions, but I couldn't find any

that drew them with proportional areas.

-Michael

--
Michael Friendly Email: friendly AT yorku DOT ca
Professor, Psychology Dept.
York University  Voice: 416 736-5115 x66249 Fax: 416 736-5814
4700 Keele StreetWeb:   http://www.datavis.ca
Toronto, ONT  M3J 1P3 CANADA

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] A question on Programming

2011-11-11 Thread Steve Lianoglou
Hi,

On Fri, Nov 11, 2011 at 7:26 AM, Oliver oli...@first.in-berlin.de wrote:
 Christofer Bogaso bogaso.christofer at gmail.com writes:
 [...]
 Here my question is, is there any speed reduction if I put them within a
 function (I think there may be some speed reduction at least within
 for-loop, because that loop needs to call that function many times),
 relative to if I used that group of codes as-it-is in many places?
 [...]


 You did not asked for it, you may know it, or may not know it:
 if you use apply functions and other vector oriented functions,
 this can bring you a huge speedup, compared to a for-loop.

I don't think this is exactly true.

I believe the majority *apply functions are really more-or-less just
sugar that boil down to having a for loop buried within them.

This is just to say that I'm not sure how much of a speedup anybody
should expect switching from a `for` loop to `*apply` function
(assuming you aren't growing an array w/in each iteration of a for
loop, for example), but if you can replace the for/*apply loop with a
few lines of vectorized functions, then yeah -- you can expect big
gains.

-- 
Steve Lianoglou
Graduate Student: Computational Systems Biology
 | Memorial Sloan-Kettering Cancer Center
 | Weill Medical College of Cornell University
Contact Info: http://cbio.mskcc.org/~lianos/contact

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] optim seems to be finding a local minimum

2011-11-11 Thread Ravi Varadhan
Hi Dimitri,

Your problem has little to do with local versus global optimum.  You can 
convince yourself that the solution you got is not even a local optimum by 
checking the gradient at the solution.

The main issue is that your objective function is not differentiable 
everywhere.  So, you have 2 options: either you use a smooth objective function 
(e.g. squared residuals) or you use an optimization algorithm than can handle 
non-smooth objective function.

Here I show that your problem is well solved by the `nmkb' function (a 
bound-constraints version of Nelder-Mead simplex method) from the dfoptim 
package.

library(dfoptim)

 myopt2 - nmkb(fn=myfunc, par=c(0.1,max(IV)), lower=0)
 myopt2
$par
[1] 8.897590e-01 9.470163e+07

$value
[1] 334.1901

$feval
[1] 204

$restarts
[1] 0

$convergence
[1] 0

$message
[1] Successful convergence

Then, there is also the issue of properly scaling your function, because it is 
poorly scaled.   Look how different the 2 parameters are - they are 7 orders of 
magnitude apart.  You are really asking for trouble here.

Hope this is helpful,
Ravi.

---
Ravi Varadhan, Ph.D.
Assistant Professor,
Division of Geriatric Medicine and Gerontology School of Medicine Johns Hopkins 
University

Ph. (410) 502-2619
email: rvarad...@jhmi.edumailto:rvarad...@jhmi.edu


[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] Predicting x from y

2011-11-11 Thread Schreiber, Stefan

Dear list members,

I have just a quick question:

I fitted a non-linear model y=a/x+b to describe my data (x=temperature and 
y=damage in %) and it works really nicely (see example below). I have 7 
different species and 8 individuals per species. I measured damage for each 
individual per species at 4 different temperatures (e.g. -5, -10, -20, -40). 
Using the individuals per species, I fitted one model per species. Now I'd like 
to use the fitted model to go back and predict the temperature that causes 50% 
damage (and it's error). Basically, it pretty easy by just rearranging the 
formula to x=a/(y-b). But that way I don't get a measure of that temperature's 
error, do I? Can I use the residual standard error that R gave me for the 
non-linear model fit? Or do I have to fit 8 lines (each individual) per 
species, calculate x based on the 8 individuals and then take the mean?

Unfortunately, dose.p from the MASS package doesn't work for non-linear models. 
When I take the log(abs(x)) the relationship becomes not satisfactory linear 
either.

Any suggestions are highly appreciated!

Thank you!
Stefan

EXAMPLE for species #1:

y.damage-c(5.7388985,1.7813519,3.7321461,2.9671031,
0.3223196,0.3207941,-1.4197658,-5.3472160,
41.1826677,29.3115243,31.3208841,35.3934115,
58.5848778,31.1541049,42.2983479,27.0615648,
64.1037728,54.7003353,66.7317044,65.4725881,
72.5755056,67.2683495,64.8717942,65.9603073,
75.0762273,56.7041960,60.0049429,70.0286506,
73.2801947,72.7015642,75.0944694,81.0361280)

x.temp-c(-5,-5,-5,-5,-5,-5,-5,-5,-10,-10,-10,-10,-10,-10,-10,
-10,-20,-20,-20,-20,-20,-20,-20,-20,-40,-40,-40,-40,-40,
-40,-40,-40)

nls(y.damage~a/x.temp+b,start=list(a=400,b=80))
plot(y.damage~x.temp,xlab='Temperature',ylab='Damage [%]')
curve(409.61/x+81.84,from=min(x.temp),to=max(x.temp),add=T)
 




[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] barplot names.arg

2011-11-11 Thread David L Carlson
That line in the documentation refers to groups (i.e. columns) on each row
of data. C.1 and C.2 are your groups and if you use the beside=TRUE option
C.1 and C.2 will appear side-by-side with a single labels on the axis
(instead of stacked one atop the other. You could rearrange your data so
that columns C.1 and C.2 in the duplicate weekdays became C3. And C.4 in the
grouped bar graph and the total number of rows in the matrix would be 4.
Using beside=TRUE would give you four bars for Mon/Tue and two bars for
Thur/Fri. But you would not have C.1/C.2 and C.3/C.4 stacked.


--
David L Carlson
Associate Professor of Anthropology
Texas AM University
College Station, TX 77843-4352


-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On
Behalf Of Diviya Smith
Sent: Thursday, November 10, 2011 11:10 PM
To: R. Michael Weylandt
Cc: r-help@r-project.org
Subject: Re: [R] barplot names.arg

The documentation for bar plots includes -

names.arga vector of names to be plotted below each bar or *group of bars*.
If this argument is omitted, then the names are taken from the names
attribute
of height if this is a vector, or the column names if it is a matrix

Thanks for your suggestion. However, I am still not sure how to group all
the entries for Mon and include one label. Your solution sort of does
that but it actually does not print the individual values . I would like to
show that there is a lot of variation in C.1 for Mon.




On Thu, Nov 10, 2011 at 11:32 PM, R. Michael Weylandt 
michael.weyla...@gmail.com wrote:

 I'm not sure you can do that data aggregation in barplot directly (a
 quick skim doesn't reveal anything in the documentation that suggests
 it to me, thought I might have missed it) though I think this does
 what you are talking about:

 barplot(sapply(unique(rownames(mdat)), function(n)
 colSums(mdat[n,,drop=F])), col = rainbow(2))

 Michael

 On Thu, Nov 10, 2011 at 9:21 PM, Diviya Smith diviya.sm...@gmail.com
 wrote:
  Hello there,
 
  I have a question regarding bar plots. I am trying to plot the data from
  the following matrix as a barplot -
 
  # input data
  mdat - matrix(c(0.1,0.9,0.9,0.1,0.5,0.5,0.45,1-0.45,0.6,0.4,0.8,0.2),
 nrow
  = 6, ncol=2, byrow=TRUE,
  +dimnames = list(c(Mon, Mon, Tues, Tues,
Thurs,
  Friday),
  +c(C.1, C.2)))
 
  # plot
 
  barplot(t(as.matrix(mdat)), col=rainbow(2), ylab=Sales, names.arg =
  rownames(mdat), border=NA, cex.names=0.5)
 
 
  I am using names.arg to print the label for the bars. However, I would
 like
  to group all the entries for the Mon and print the label only once. Is
  there a way to do this? The documentation for barplots suggests that
this
  can be done but I was not able to figure it out. Please help.
 
 
  Thanks in advance,
 
  Diviya
 
 [[alternative HTML version deleted]]
 
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
 


[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Upgrade R?

2011-11-11 Thread Uwe Ligges



On 10.11.2011 17:53, Kevin Burton wrote:

The problem with this documentation is two-fold. One it seems to concentrate
on building from source which I don't need. Two it doesn't address the
upgade. I have a number of packages and so I need to do what has been
suggested and install the latest version *first*. Then copy the libraries
(packages). Then uninstall the previous version. It is on this last step
that I am stuck on right now. The last link on uninstalling R manually was
what I needed. Thank you.


Just delete the directory.

Uwe Ligges








Kevin



From: jose Bartolomei [mailto:surfpr...@hotmail.com]
Sent: Thursday, November 10, 2011 10:19 AM
To: rkevinbur...@charter.net; R Help
Subject: RE: [R] Upgrade R?



Hi,
Don't know if this will help you but...
In my short experience and following the guidelines you should first
uninstall R.

http://cran.r-project.org/doc/manuals/R-admin.html#Installing-R-under-Window
s

Unistall it from the Windows control panel.

The old R version libraries file will be kept on machine.
For example : C:\Program Files\R\R-2.13.0\library

Then install the new version via:
http://cran.r-project.org/doc/manuals/R-admin.html#Installing-R-under-Window
s

You can copy/paste libraries from the old version R library file to the new
one.
C:\Program Files\R\R-2.14.0\library

There is too an function named:
?update.packagee

If above was what you did, then there is a post on Uninstalling R manually:

http://learnserver.csd.univie.ac.at/rcomwiki/doku.php?id=wiki:uninstalling_r
_manually

Regards,
Jose



From: rkevinbur...@charter.net
To: r-help@r-project.org
Date: Thu, 10 Nov 2011 09:07:20 -0600
Subject: Re: [R] Upgrade R?

Since apparently there is no one familiar with this error message let me
rephrase the question. Is there a 'manual' process to fully remove a

version

of 'R' from my machine? This is a Window PC running Windows 7.



Thank you.



Kevin



From: Kevin Burton [mailto:rkevinbur...@charter.net]
Sent: Monday, November 07, 2011 2:23 PM
To: 'r-help@r-project.org'
Subject: Upgrade R?



I am trying to upgrade to R 2.14 from R 2.13.1 I have compied all the
libraries from the 'library' directory in my existing installation

(2.13.1)

to the installed R 2.14. Now I want to uninstall the old installation (R
2.13.1) and I get the error:



Internal Error: Cannot find utCompiledCode record for this version of the
uninstaller.



Any ideas?



Kevin


[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide

http://www.R-project.org/posting-guide.html

and provide commented, minimal, self-contained, reproducible code.



[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] optim seems to be finding a local minimum

2011-11-11 Thread Dimitri Liakhovitski
Thank you very much to everyone who replied!
As I mentioned - I am not a mathematician, so sorry for stupid
comments/questions.
I intuitively understand what you mean by scaling. While the solution
space for the first parameter (.alpha) is relatively compact (probably
between 0 and 2), the second one (.beta) is all over the place -
because it is a function of IV (input vector). And that's, probably,
my main challenge - that I am trying to write a routine for different
possible IVs that I might be facing (they may be in hundreds, in
thousands, in millions). Should I be rescaling the IV somehow (e.g.,
by dividing it by its max) - or should I do something with the
parameter .beta inside my function?

So far, I've written a loop over many different starting points for
both parameters. Then, I take the betas around the best solution so
far, split it into smaller steps for beta (as starting points) and
optimize again for those starting points. What disappoints me is that
even when I found a decent solution (the minimized value of 336) it
was still worse than the Solver solution!

And I am trying to prove to everyone here that we should do R, not Excel :-)

Thanks again for your help, guys!
Dimitri


On Fri, Nov 11, 2011 at 9:10 AM, John C Nash nas...@uottawa.ca wrote:
 I won't requote all the other msgs, but the latest (and possibly a bit 
 glitchy) version of
 optimx on R-forge

 1) finds that some methods wander into domains where the user function fails 
 try() (new
 optimx runs try() around all function calls). This includes L-BFGS-B

 2) reports that the scaling is such that you really might not expect to get a 
 good solution

 then

 3) Actually gets a better result than the

 xlf-myfunc(c(0.888452533990788,94812732.0897449))
 xlf
 [1] 334.607


 with Kelley's variant of Nelder Mead (from dfoptim package), with

 myoptx
  method                        par       fvalues fns  grs itns conv  KKT1
 4 LBFGSB                     NA, NA 8.988466e+307  NA NULL NULL     NA
 2 Rvmmin           0.1, 200186870.6      25593.83  20    1 NULL    0 FALSE
 3 bobyqa 6.987875e-01, 2.001869e+08      1933.229  44   NA NULL    0 FALSE
 1   nmkb 8.897590e-01, 9.470163e+07      334.1901 204   NA NULL    0 FALSE
   KKT2 xtimes  meths
 4    NA   0.01 LBFGSB
 2 FALSE   0.11 Rvmmin
 3 FALSE   0.24 bobyqa
 1 FALSE   1.08   nmkb

 But do note the terrible scaling. Hardly surprising that this function does 
 not work. I'll
 have to delve deeper to see what the scaling setup should be because of the 
 nature of the
 function setup involving some of the data. (optimx includes parscale on all 
 methods).

 However, original poster DID include code, so it was easy to do a quick 
 check. Good for him.

 JN

 ## Comparing this solution to Excel Solver solution:
 myfunc(c(0.888452533990788,94812732.0897449))

 -- Dimitri Liakhovitski marketfusionanalytics.com

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.




-- 
Dimitri Liakhovitski
marketfusionanalytics.com

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] optim seems to be finding a local minimum

2011-11-11 Thread John C Nash
Some tips:

1) Excel did not, as far as I can determine, find a solution. No point seems to 
satisfy
the KKT conditions (there is a function kktc in optfntools on R-forge project 
optimizer.
It is called by optimx).

2) Scaling of the input vector is a good idea given the seeming wide range of 
values. That
is, assuming this can be done. If the function depends on the relative values 
in the input
vector rather than magnitude, this may explain the trouble with your function. 
That is, if
the function depends on the relative change in the input vector and not its 
scale, then
optimizers will have a lot of trouble if the scale factor for this vector is 
implicitly
one of the optimization parameters.

3) If you can get the gradient function you will almost certainly be able to do 
better,
especially in finding whether you have a minimum i.e., null gradient, positive 
definite
Hessian. When you have gradient function, kktc uses Jacobian(gradient) to get 
the Hessian,
avoiding one level of digit cancellation.

JN


On 11/11/2011 10:20 AM, Dimitri Liakhovitski wrote:
 Thank you very much to everyone who replied!
 As I mentioned - I am not a mathematician, so sorry for stupid
 comments/questions.
 I intuitively understand what you mean by scaling. While the solution
 space for the first parameter (.alpha) is relatively compact (probably
 between 0 and 2), the second one (.beta) is all over the place -
 because it is a function of IV (input vector). And that's, probably,
 my main challenge - that I am trying to write a routine for different
 possible IVs that I might be facing (they may be in hundreds, in
 thousands, in millions). Should I be rescaling the IV somehow (e.g.,
 by dividing it by its max) - or should I do something with the
 parameter .beta inside my function?
 
 So far, I've written a loop over many different starting points for
 both parameters. Then, I take the betas around the best solution so
 far, split it into smaller steps for beta (as starting points) and
 optimize again for those starting points. What disappoints me is that
 even when I found a decent solution (the minimized value of 336) it
 was still worse than the Solver solution!
 
 And I am trying to prove to everyone here that we should do R, not Excel :-)
 
 Thanks again for your help, guys!
 Dimitri
 
 
 On Fri, Nov 11, 2011 at 9:10 AM, John C Nash nas...@uottawa.ca wrote:
 I won't requote all the other msgs, but the latest (and possibly a bit 
 glitchy) version of
 optimx on R-forge

 1) finds that some methods wander into domains where the user function fails 
 try() (new
 optimx runs try() around all function calls). This includes L-BFGS-B

 2) reports that the scaling is such that you really might not expect to get 
 a good solution

 then

 3) Actually gets a better result than the

 xlf-myfunc(c(0.888452533990788,94812732.0897449))
 xlf
 [1] 334.607


 with Kelley's variant of Nelder Mead (from dfoptim package), with

 myoptx
  methodpar   fvalues fns  grs itns conv  KKT1
 4 LBFGSB NA, NA 8.988466e+307  NA NULL NULL NA
 2 Rvmmin   0.1, 200186870.6  25593.83  201 NULL0 FALSE
 3 bobyqa 6.987875e-01, 2.001869e+08  1933.229  44   NA NULL0 FALSE
 1   nmkb 8.897590e-01, 9.470163e+07  334.1901 204   NA NULL0 FALSE
   KKT2 xtimes  meths
 4NA   0.01 LBFGSB
 2 FALSE   0.11 Rvmmin
 3 FALSE   0.24 bobyqa
 1 FALSE   1.08   nmkb

 But do note the terrible scaling. Hardly surprising that this function does 
 not work. I'll
 have to delve deeper to see what the scaling setup should be because of the 
 nature of the
 function setup involving some of the data. (optimx includes parscale on all 
 methods).

 However, original poster DID include code, so it was easy to do a quick 
 check. Good for him.

 JN

 ## Comparing this solution to Excel Solver solution:
 myfunc(c(0.888452533990788,94812732.0897449))

 -- Dimitri Liakhovitski marketfusionanalytics.com

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

 
 


__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] Fwd: Use of R for VECM

2011-11-11 Thread vramaiah

- Forwarded Message -
From: vrama...@neo.tamu.edu
To: bernhard pfaff bernhard.pf...@pfaffikus.de
Sent: Friday, November 11, 2011 9:03:11 AM GMT -06:00 US/Canada Central
Subject: Use of R for VECM

Hello Fellow R'ers
I am a new user of R and I am applying it for solving Bi-Variate (Consumption 
and Output) VECM with Co-Integration (I(1)) with three lags on Consumption and 
Output expressed aa lags of differences.  R Code and one segment fo the output 
(other parts of the output are repeatitive) are as follows.
 us=read.table(usdata.1995-2009.txt,header=T)
 sjd-us[,c(Y,C)]
 sjd.vecm1 - ca.jo(sjd, ecdet='const', type=eigen, K=3, spec=longrun,
+ season=4)
 sjd.vecm2 - ca.jo(sjd, ecdet='const', type=eigen, K=3, spec=transitory,
+ season=4)
 sjd.vecm.ols1 - cajools(sjd.vecm1)
 sjd.vecm.ols2 - cajools(sjd.vecm2)
 summary(sjd.vecm.ols1)
Response Y.d :

Call:
lm(formula = substitute(Y.d), data = data.mat)

Residuals:
   Min 1Q Median 3QMax 
-0.0049787 -0.0012948  0.703  0.0009653  0.0063192 

Coefficients:
   Estimate Std. Error t value Pr(|t|)
sd1  -0.0002012  0.0007653  -0.263 0.793820
sd2   0.0013339  0.0007616   1.752 0.086378 .  
sd3   0.0007372  0.0007947   0.928 0.358348
Y.dl1-0.2215246  0.1600532  -1.384 0.172875
C.dl1 0.9220846  0.1646314   5.601 1.08e-06 ***
Y.dl2-0.1600219  0.1273838  -1.256 0.215245
C.dl2 0.4712112  0.2124175   2.218 0.031401 *  
Y.l3 -0.3227708  0.0881079  -3.663 0.000631 ***
C.l3  0.2376579  0.0688854   3.450 0.001194 ** 
constant  0.4707624  0.1182284   3.982 0.000236 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.002013 on 47 degrees of freedom
Multiple R-squared: 0.7053, Adjusted R-squared: 0.6426 
F-statistic: 11.25 on 10 and 47 DF,  p-value: 1.707e-09 



I am trying to decipher the output and have the following questions.  
1.  If you have published a book with explanations of the output, I would like 
to acquire one.  Please advise me where I can get one.  It might me and you lot 
of time.  In the meantime
2.  Are Y.dl1, C.dl1, Y.dl2, C.dl2, Y.l3 and C.l3 are the coefficients I am 
looking for?
3.  My equition requires Y.l1 and C.l1.  Is there a way I could change from 
Y.l3 and C.l3 
4.  I need to use these coefficients for forecasting.  What command I use to 
extract the coefficients from this program?
5.  What are sd1, sd2, sd3?  I did not input any seasonal dummies and yet the 
output seems to have included them.  Does it affect the value of the 
coefficients?  Is there a better way to reframe the command to avoid them?
I hope it is not too much to ask for your help.
I thank you in advance.
Ram Ramaiah

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] performance of adaptIntegrate vs. integrate

2011-11-11 Thread Ravi Varadhan
The integrand is highly peaked.  It is approximately an impulse function where 
much of the mass is concentrated at a very small interval.  Plot the function 
and see for yourself.  This is the likely cause of the problem.

Other types of integrands where you could experience problems are: integrands 
with singularity at either limit and slowly decaying oscillatory integrands.  
As to why integrate performs better than adaptIntegrate in this situation, I 
don't know.  You have to study the two implementations.   Wynn's epsilon 
algorithm is an extrapolation method for improving the convergence of a 
sequence.  This could be an explanation for the better performance, but I 
cannot say for sure.

Hope this is helpful,
Ravi
---
Ravi Varadhan, Ph.D.
Assistant Professor,
Division of Geriatric Medicine and Gerontology School of Medicine Johns Hopkins 
University

Ph. (410) 502-2619
email: rvarad...@jhmi.edumailto:rvarad...@jhmi.edu


[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] drawing ellipses in R

2011-11-11 Thread Greg Snow
Those formulas are the standard way to convert from polar coordinates to 
Euclidean coordinates.  The polar coordinates are 'r' which is the radius or 
distance from the center point and 'theta' which is the angle (0 is pointing in 
the positive x direction).

If r is constant and theta coveres a full cycle then you will get a circle.

-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
Behalf Of mms...@comcast.net
Sent: Monday, October 31, 2011 10:50 PM
To: r-h...@stat.math.ethz.ch
Subject: [R] drawing ellipses in R

Hello, 

I have been following the thread dated Monday, October 9, 2006 when Kamila 
Naxerova asked a question about plotting elliptical shapes. Can you explain the 
equations for X and Y. I believe they used the parametric form of x and y (x=r 
cos(theta), y=r sin(theta). I don't know what r is here ? Can you explain 1)the 
origin of these equations and 2) what is r? 

Sincerely, 
Mary A. Marion 

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] multivariate random variable

2011-11-11 Thread Greg Snow
Your question is a bit too general to give a useful answer.  One possible 
answer to your question is:

 mrv - matrix( runif(1000), ncol=10 )

Which generates multivariate random observations, but is unlikely to be what 
you are really trying to accomplish.  There are many tools for generating 
multivariate random data including Metropolis-Hastigns, Gibbs sampling, 
rejection sampling, conditional generation, copulas, and many others, which one 
will be best (or which combination will be best) depends on what you are 
actually trying to accomplish.

-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
Behalf Of Anera Salucci
Sent: Tuesday, November 01, 2011 5:22 AM
To: r-help@r-project.org
Subject: [R] multivariate random variable

Dear All,
 
How can I generate multivariate random variable (not multivariate normal )
 
I am in urgent
[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Building a statically-linked extension?

2011-11-11 Thread Prof Brian Ripley

On Fri, 11 Nov 2011, Tyler Pirtle wrote:




On Fri, Nov 11, 2011 at 2:09 AM, Prof Brian Ripley rip...@stats.ox.ac.uk
wrote:
  On Thu, 10 Nov 2011, Tyler Pirtle wrote:

Hi there,

I'm writing an R extension that has a C component
that relies on two third
party libraries that I'm bundling
with the extension.

I'd like to statically link my resulting extension
so that I can rely on
the bundled versions of the libraries I'm
distributing, so there's two questions -

1, does R allow statically linked C extensions to be
used at runtime?


Yes

  2, are there any standard ways of having R build my
  extension statically?


Yes.

If you want to know more, follow the posting guide: this is not an
R-help topic, the missing 'at a minimum' information is vital 


Apologies for spamming R-help, I suppose this should have gone to R-devel. I
can't ever remember the difference between
the two. So, Professor, if I re-post to R-devel would you be willing to
reveal any more detail for me? ;)

I also apologize for a lack of additional information, I'm happy to provide
you with any, but I'm not exactly sure
what else I can offer..I just want to build a statically linked R extension,
and it doesn't seem to be covered in Writing R
Extensions. You seem to know the answer, any pointers would be greatly
appreciated.


You will need to tell us your OS, for a start.




Tyler

 

Thanks,


Tyler

       [[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained,
reproducible code.


  --
  Brian D. Ripley,                  rip...@stats.ox.ac.uk
  Professor of Applied Statistics,
   http://www.stats.ox.ac.uk/~ripley/
  University of Oxford,             Tel:  +44 1865 272861 (self)
  1 South Parks Road,                     +44 1865 272866 (PA)
  Oxford OX1 3TG, UK                Fax:  +44 1865 272595






--
Brian D. Ripley,  rip...@stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Building a statically-linked extension?

2011-11-11 Thread Tyler Pirtle
On Fri, Nov 11, 2011 at 2:09 AM, Prof Brian Ripley rip...@stats.ox.ac.ukwrote:

 On Thu, 10 Nov 2011, Tyler Pirtle wrote:

  Hi there,

 I'm writing an R extension that has a C component that relies on two third
 party libraries that I'm bundling
 with the extension.

 I'd like to statically link my resulting extension so that I can rely on
 the bundled versions of the libraries I'm
 distributing, so there's two questions -

 1, does R allow statically linked C extensions to be used at runtime?


 Yes


  2, are there any standard ways of having R build my extension statically?


 Yes.

 If you want to know more, follow the posting guide: this is not an R-help
 topic, the missing 'at a minimum' information is vital 


Apologies for spamming R-help, I suppose this should have gone to R-devel.
I can't ever remember the difference between
the two. So, Professor, if I re-post to R-devel would you be willing to
reveal any more detail for me? ;)

I also apologize for a lack of additional information, I'm happy to provide
you with any, but I'm not exactly sure
what else I can offer..I just want to build a statically linked R
extension, and it doesn't seem to be covered in Writing R
Extensions. You seem to know the answer, any pointers would be greatly
appreciated.


Tyler




 Thanks,


 Tyler

[[alternative HTML version deleted]]

 __**
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/**listinfo/r-helphttps://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/**
 posting-guide.html http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


 --
 Brian D. Ripley,  rip...@stats.ox.ac.uk
 Professor of Applied Statistics,  
 http://www.stats.ox.ac.uk/~**ripley/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


[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] Jordan Form of a matrix

2011-11-11 Thread Arnau Mir Torres
Hello.

Is it possible to find the Jordan Form of a matrix with R?


Arnau.

Arnau Mir Torres
Edifici A. Turmeda
Campus UIB
Ctra. Valldemossa, km. 7,5
07122 Palma de Mca.
tel: (+34) 971172987
fax: (+34) 971173003
email: arnau@uib.es
URL: http://dmi.uib.es/~arnau

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Fwd: Use of R for VECM

2011-11-11 Thread R. Michael Weylandt michael.weyla...@gmail.com
I suspect that part of the hesitancy in replying to your query is that your 
academic e-mail and subsequent google-ability suggest you are a student in the 
economics department and that this might be homework. 

In the meanwhile, you might want to look into Pfaff's book on cointegration and 
time series analysis in R. 

Michael

On Nov 11, 2011, at 10:06 AM, vrama...@neo.tamu.edu wrote:

 
 - Forwarded Message -
 From: vrama...@neo.tamu.edu
 To: bernhard pfaff bernhard.pf...@pfaffikus.de
 Sent: Friday, November 11, 2011 9:03:11 AM GMT -06:00 US/Canada Central
 Subject: Use of R for VECM
 
 Hello Fellow R'ers
 I am a new user of R and I am applying it for solving Bi-Variate (Consumption 
 and Output) VECM with Co-Integration (I(1)) with three lags on Consumption 
 and Output expressed aa lags of differences.  R Code and one segment fo the 
 output (other parts of the output are repeatitive) are as follows.
 us=read.table(usdata.1995-2009.txt,header=T)
 sjd-us[,c(Y,C)]
 sjd.vecm1 - ca.jo(sjd, ecdet='const', type=eigen, K=3, spec=longrun,
 + season=4)
 sjd.vecm2 - ca.jo(sjd, ecdet='const', type=eigen, K=3, spec=transitory,
 + season=4)
 sjd.vecm.ols1 - cajools(sjd.vecm1)
 sjd.vecm.ols2 - cajools(sjd.vecm2)
 summary(sjd.vecm.ols1)
 Response Y.d :
 
 Call:
 lm(formula = substitute(Y.d), data = data.mat)
 
 Residuals:
   Min 1Q Median 3QMax 
 -0.0049787 -0.0012948  0.703  0.0009653  0.0063192 
 
 Coefficients:
   Estimate Std. Error t value Pr(|t|)
 sd1  -0.0002012  0.0007653  -0.263 0.793820
 sd2   0.0013339  0.0007616   1.752 0.086378 .  
 sd3   0.0007372  0.0007947   0.928 0.358348
 Y.dl1-0.2215246  0.1600532  -1.384 0.172875
 C.dl1 0.9220846  0.1646314   5.601 1.08e-06 ***
 Y.dl2-0.1600219  0.1273838  -1.256 0.215245
 C.dl2 0.4712112  0.2124175   2.218 0.031401 *  
 Y.l3 -0.3227708  0.0881079  -3.663 0.000631 ***
 C.l3  0.2376579  0.0688854   3.450 0.001194 ** 
 constant  0.4707624  0.1182284   3.982 0.000236 ***
 ---
 Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
 
 Residual standard error: 0.002013 on 47 degrees of freedom
 Multiple R-squared: 0.7053,Adjusted R-squared: 0.6426 
 F-statistic: 11.25 on 10 and 47 DF,  p-value: 1.707e-09 
 
 
 
 I am trying to decipher the output and have the following questions.  
 1.  If you have published a book with explanations of the output, I would 
 like to acquire one.  Please advise me where I can get one.  It might me and 
 you lot of time.  In the meantime
 2.  Are Y.dl1, C.dl1, Y.dl2, C.dl2, Y.l3 and C.l3 are the coefficients I am 
 looking for?
 3.  My equition requires Y.l1 and C.l1.  Is there a way I could change from 
 Y.l3 and C.l3 
 4.  I need to use these coefficients for forecasting.  What command I use to 
 extract the coefficients from this program?
 5.  What are sd1, sd2, sd3?  I did not input any seasonal dummies and yet 
 the output seems to have included them.  Does it affect the value of the 
 coefficients?  Is there a better way to reframe the command to avoid them?
 I hope it is not too much to ask for your help.
 I thank you in advance.
 Ram Ramaiah
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] rbind.data.frame drops attributes for factor variables

2011-11-11 Thread Jeff Newmiller
As the doctor says, if it hurts don't do that.

A factor is a sequence of integers with a corresponding list of character 
strings. Factors in two separate vectors can and usually do map the same 
integer to different strings, and R cannot tell how you want that resolved.

Convert these columns to character before combining them, and only convert to 
factor when you have all of your possibilities present (or you specify them in 
the creation of the factor vector).
---
Jeff NewmillerThe .   .  Go Live...
DCN:jdnew...@dcn.davis.ca.usBasics: ##.#.   ##.#.  Live Go...
  Live:   OO#.. Dead: OO#..  Playing
Research Engineer (Solar/BatteriesO.O#.   #.O#.  with
/Software/Embedded Controllers)   .OO#.   .OO#.  rocks...1k
--- 
Sent from my phone. Please excuse my brevity.

Sammy Zee szee2...@gmail.com wrote:

Hi all,

When I use rbind() or rbind.data.frame() to add a row to an existing
dataframe, it appears that attributes for the column of type factor
are
dropped. I see the following post with same problem. However i did not
see
any reply to the following posting offering a solution. Could someone
please help.

http://r.789695.n4.nabble.com/rbind-data-frame-drops-attributes-for-factor-variables-td919575.html

Thanks,
Sanjay

   [[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] mc.cores and computer settings on osx and linux

2011-11-11 Thread ivo welch
for the googleable r-help archives, I thought I would post what I
wrote into my .Rprofile to automatically set some system information.
the most relevant aspect is the determination of mc.cores.  this is
useful when users want to use the parallel package

  options(uname= system(uname, intern=TRUE))

  options(os= if (getOption(uname)==Darwin) osx else linux)

  if ((getOption(os) != osx)  (getOption(os) != linux))
stop(You need to set options yourself.  I only grok linux and osx\n)

  options(mc.cores= as.numeric(if (getOption(os)==osx)
system(sysctl -n hw.ncpu, intern=TRUE) else system(grep 'core id'
/proc/cpuinfo | sort | uniq | wc -l, intern=TRUE)))

  options(hostname= system(hostname, intern=TRUE))

[below, I am also posting my current wrapper for the parallel library.
 I know it is amateurish, but it may be useful for novices exploring
parallel calculations.  it is a friendlier face.]

[gripes: R is powerful, and the team that maintains it are saints.
But R is not friendly.  it lacks the ability to turn off recycling for
enhanced error detection.  it does not throw clear errors when one
accesses a non-existing column in a data frame.  it does not print out
the user program line number where an error occurred.  it lacks an
end-user documentation system [like POD], though it does have good
package documentation.  it does have some unexpected behavior:
mymatrix[1:2,] is a matrix, but mymatrix[1:1,] is a numeric.  huh?
data.table is necessary for reasonably fast data manipulation, but
data.table giveth and taketh.  it has some really strange unexpected
behavior---mydatatable[,1] is not the second column, as one would
expect it to be.  yes, it is documented, but syntax should be as
expected.]

/iaw


Ivo Welch (ivo.we...@gmail.com)




###
### these R functions are very type-limited wrappers for
### by()-like operations, using the multicore library.  this
### means effort-less multi-CPU calculations.
###
### the user functions MUST return a numeric scalar, a vector, a matrix, or
### a data frame.  to enhance speed, internally the user function is
### wrapped, too.
###
### the output is ONE matrix, whose row-names are the categories.


check.output - function( mc.rv ) {
  ## check that we have a list of matrices, and that each matrix has
the same number of columns
  numofcols - (-1)
  for (i in 1:length(mc.rv)) {
if (is.null(mc.rv[[i]])) next;
if (! (is.matrix(mc.rv[[i]])|is.data.frame(mc.rv[[i]])) )
  abort(iaw-mc.R:check.output: Element, i, is not a
matrix/dataframe, but a , whatis(mc.rv[[i]]))
if (numofcols0) numofcols - ncol(mc.rv[[i]])
if (numofcols0) next
if (ncol(mc.rv[[i]]) != numofcols) {
  print(head(mc.rv[[i]]))
  abort(iaw-mc.R:check.output: Element, i, should have,
numofcols, columns, but has, ncol(mc.rv[[i]]), columns instead.)
}
  }
}

add.by.names - function( mc.rv ) {
  for (i in 1:length(mc.rv))
if (!is.null(mc.rv[[i]])) row.names(mc.rv[[i]]) - rep(
names(mc.rv)[i], nrow(mc.rv[[i]]) )
  mc.rv
}


.mc.by - function(lcapplyversion, data, INDICES, FUN, ...) {
  si - split(1:nrow(data), INDICES)

  ## input = set of row indexes ; output = one row in a matrix or data
frame, that can be stacked up
  FUN.ON.ROWS - function(.index, ...) { rv - FUN(data[.index,],
...); if (is.null(rv)) rv else if (is.vector(rv)) matrix(rv, nrow=1)
else rv }
  soln - lcapplyversion( si, FUN.ON.ROWS, ... )
  check.output(soln)

  rv - do.call(rbind, add.by.names(soln))
  if (is.null(rv)) { print(head(soln)); abort(Sorry, but in .mc.by,
the rv is null!\n) }
  if (ncol(rv)==1) {
nm - rownames(rv)
rv - as.vector(rv)
names(rv) - nm
  }
  rv
}

mc.by - function(data, INDICES, FUN, ...) { .mc.by(mclapply, data,
INDICES, FUN, ...) }
oc.by - function(data, INDICES, FUN, ...) { .mc.by(lapply, data,
INDICES, FUN, ...) }


mc.byallrows - function(data, FUN, ...) {
  si - as.list(1:nrow(data))  ## a little faster than the split for
large data sets
  FUN.ON.ROWS - function(.index, ...) { rv - FUN(data[.index,],
...); if (is.null(rv)) rv else if (is.vector(rv)) matrix(rv, nrow=1)
else rv }

  soln - mclapply( si, FUN.ON.ROWS, ..., mc.cores= 4 )
  check.output(soln)
  rv - do.call(rbind, soln)  ## omits naming.
  if (ncol(rv)==1) rv - as.vector(rv)
  rv
}


if (0) {
  function.sample - function(d) cbind(d$x+d$y, d$x, d$y)
  function.sample.simpler - function(d) (d$x+d$y)

  d - data.frame( i=c( rep(1,2), rep(2,3), rep(3,4) ), x=rnorm(9), y=rnorm(9) )

  report - function( text2print, f.output ) {
cat(\n\n, text2print, :\n); print(f.output); cat(\n\n)
  }

  report( the original R by() function, by( d, d$i, function.sample ))
  report( wrappled multicore by mc.by with user function returning
scalar, mc.by( d, d$i, function.sample.simpler ))
  report( wrappled multicore by mc.by with user function returning
vector, mc.by( d, d$i, function.sample 

[R] Random-walk Metropolis-Hasting

2011-11-11 Thread Gyanendra Pokharel
Following is my code, can some one help on the error at the bottom?

 mh-function(iterations,alpha,beta){
+ data-read.table(epidemic.txt,header = TRUE)
+ attach(data, warn.conflicts = F)
+ k-97
+ d - (sqrt((x-x[k])^2 + (y-y[k])^2))
+ p - 1-exp(-alpha*d^(-beta))
+ p.alpha-1 - exp(-3*d^(-beta))
+ p.beta - 1 - exp(alpha*d^(-2))
+ iterations-1000
+ mu.lambda - c(0,0);s.lambda - c(100,100)
+ prop.s - c(0.1,0.1)
+ lambda - matrix(nrow=iterations, ncol=2)
+ acc.prob-0
+ current.lambda - c(0,0)
+ for(t in 1:iterations){
+ prop.lambda - rnorm(2,current.lambda,prop.s)
+ a - p.beta/p.alpha
*(dnorm(prop.lambda,mu.lambda,s.lambda))*dnorm(current.lambda,mu.lambda,s.lambda)
+ accept - min(1,a)
+ u-runif(1)
+ if(u[t]=accept[t]){
+ current.lambda - prop.lambda
+ acc.prob - acc.prob +1
+ }
+ lambda[t,] - current.lambda
+ }
+ lambda
+ }
 mh(1000,0,0)
Error in if (u[t] = accept[t]) { : missing value where TRUE/FALSE needed

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] Formula variable help

2011-11-11 Thread rkevinburton

I have an R script with the following applicable lines:

xshort - window(s, start=st, end=ed)
. . .
xshort - ts(xshort, frequency=1, start=1)
. . .
m1 - m2 - m3 - m4 - m5 - m6 - NULL
m1 - tslm(xshort ~ trend)

I get an error:

Error in get(dataname) : object 'xshort' not found

When I do traceback() I get:

3: get(dataname)
2: tslm(xshort ~ trend) at #19
1: model.cross.validation(l[[MEN]]$series)

Which points to the call to tslm above.  Since I am not supply 'data' to 
the tslm call (in the forecast package),, I am assuming that the code is 
dying here (in tslm):

 if (missing(data)) {
 dataname - as.character(formula)[2]
 x - get(dataname)
 data - data.frame(x)
 colnames(data) - dataname
 }

Any ideas what is failing?

Kevin

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Formula variable help

2011-11-11 Thread Joshua Wiley
Seems like it could work---can you save your script as a .txt file and
send it as an attachment, upload it online, or make a reproducible
example?

Cheers,

Josh

On Fri, Nov 11, 2011 at 9:52 AM,  rkevinbur...@charter.net wrote:

 I have an R script with the following applicable lines:

                xshort - window(s, start=st, end=ed)
 . . .
                xshort - ts(xshort, frequency=1, start=1)
 . . .
                m1 - m2 - m3 - m4 - m5 - m6 - NULL
                m1 - tslm(xshort ~ trend)

 I get an error:

 Error in get(dataname) : object 'xshort' not found

 When I do traceback() I get:

 3: get(dataname)
 2: tslm(xshort ~ trend) at #19
 1: model.cross.validation(l[[MEN]]$series)

 Which points to the call to tslm above.  Since I am not supply 'data' to
 the tslm call (in the forecast package),, I am assuming that the code is
 dying here (in tslm):

     if (missing(data)) {
         dataname - as.character(formula)[2]
         x - get(dataname)
         data - data.frame(x)
         colnames(data) - dataname
     }

 Any ideas what is failing?

 Kevin

        [[alternative HTML version deleted]]

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.




-- 
Joshua Wiley
Ph.D. Student, Health Psychology
Programmer Analyst II, ATS Statistical Consulting Group
University of California, Los Angeles
https://joshuawiley.com/

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Random-walk Metropolis-Hasting

2011-11-11 Thread Uwe Ligges



On 11.11.2011 18:48, Gyanendra Pokharel wrote:

Following is my code, can some one help on the error at the bottom?


mh-function(iterations,alpha,beta){

+ data-read.table(epidemic.txt,header = TRUE)
+ attach(data, warn.conflicts = F)
+ k-97
+ d- (sqrt((x-x[k])^2 + (y-y[k])^2))
+ p- 1-exp(-alpha*d^(-beta))
+ p.alpha-1 - exp(-3*d^(-beta))
+ p.beta- 1 - exp(alpha*d^(-2))
+ iterations-1000
+ mu.lambda- c(0,0);s.lambda- c(100,100)
+ prop.s- c(0.1,0.1)
+ lambda- matrix(nrow=iterations, ncol=2)
+ acc.prob-0
+ current.lambda- c(0,0)
+ for(t in 1:iterations){
+ prop.lambda- rnorm(2,current.lambda,prop.s)
+ a- p.beta/p.alpha
*(dnorm(prop.lambda,mu.lambda,s.lambda))*dnorm(current.lambda,mu.lambda,s.lambda)
+ accept- min(1,a)
+ u-runif(1)
+ if(u[t]=accept[t]){
+ current.lambda- prop.lambda
+ acc.prob- acc.prob +1
+ }
+ lambda[t,]- current.lambda
+ }
+ lambda
+ }

mh(1000,0,0)

Error in if (u[t]= accept[t]) { : missing value where TRUE/FALSE needed


As the error message tells you:

accept[t] is NA t  1 since

accept- min(1,a)

Uwe Ligges






[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] symbols and legend: how to harmonize point size ?

2011-11-11 Thread Patrick Giraudoux

Hi,

I was wondering if it is possible to harmonize the ouput of symbols() 
and legend() both from the graphics package.


Let us take this example:

x-runif(10)
y-runif(10)
z-runif(10)

leg-round(seq(min(z),max(z),l=4),2) # 4 values rounded up to 2 decimals 
for the legend


symbols(x,y,circles=z,inches=0.2)

legend(topright,legend=leg,pch=1,pt.cex=leg/max(leg)*2) # multiplied 
by 2 arbitrarily just to make it visible


Actually, what I want to do is to pass to pt.cex a value which would 
make the biggest circle in the legend (leg/max(leg) = 1) exactly the 
same size as the one specified in symbols (here 0.2 inches). I suppose 
this is possible using par(cin) but I cannot figure out how to do it 
properly.


Any hint appreciated,

Best,

Patrick

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] Encryption package of R

2011-11-11 Thread Jun Ji
Dear all:

Is there any encryption package of R? For instance, 3DES or AES algorithm
package or whatever...
Thanks a lot in advance.

Best regards,
Jason

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] R Startup options (MDI and SDI)

2011-11-11 Thread Gene Leynes
Here's my setup:

   - I'm on a Windows machine (I don't have full admin rights)
   - I have a folder with an *.RData file and an .RProfile file
   - I want the user to be able to start R by double clicking on the
   *.RData file

Can I specify the application start up options (like  --no-save --mdi) in
the local RProfile file?

Thanks!

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Jordan Form of a matrix

2011-11-11 Thread Mehmet Suzen
Hi Arnau,

Not aware of direct implementation. It was discussed in octave project as well 

(http://octave.1599824.n4.nabble.com/Jordan-canonical-form-td2216965.html)
It is numerically ill-conditioned to compute that.

See this http://dl.acm.org/citation.cfm?id=355912

Best,

Mehmet Süzen, PhD

Mango Solutions
data analysis that delivers
t:   +44 (0) 1249 767700
m:   +44 (0) 7517 639318
f:   +44 (0) 1249 767707
e: msu...@mango-solutions.com
w: www.mango-solutions.com


-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org]
On Behalf Of Arnau Mir Torres
Sent: 11 November 2011 17:19
To: r-help@r-project.org
Subject: [R] Jordan Form of a matrix

Hello.

Is it possible to find the Jordan Form of a matrix with R?


Arnau.

Arnau Mir Torres
Edifici A. Turmeda
Campus UIB
Ctra. Valldemossa, km. 7,5
07122 Palma de Mca.
tel: (+34) 971172987
fax: (+34) 971173003
email: arnau@uib.es
URL: http://dmi.uib.es/~arnau

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-
guide.html
and provide commented, minimal, self-contained, reproducible code.
LEGAL NOTICE
This message is intended for the use o...{{dropped:10}}

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Formula variable help

2011-11-11 Thread Kevin Burton
It seems that there is a bug in the forecast::tslm function. I have forwarded 
what I think is the bug (the call to get() should supply the argument 
'envir=parent.frame()').

Thank you.

On Nov 11, 2011, at 12:11 PM, Joshua Wiley jwiley.ps...@gmail.com wrote:

 Seems like it could work---can you save your script as a .txt file and
 send it as an attachment, upload it online, or make a reproducible
 example?
 
 Cheers,
 
 Josh
 
 On Fri, Nov 11, 2011 at 9:52 AM,  rkevinbur...@charter.net wrote:
 
 I have an R script with the following applicable lines:
 
xshort - window(s, start=st, end=ed)
 . . .
xshort - ts(xshort, frequency=1, start=1)
 . . .
m1 - m2 - m3 - m4 - m5 - m6 - NULL
m1 - tslm(xshort ~ trend)
 
 I get an error:
 
 Error in get(dataname) : object 'xshort' not found
 
 When I do traceback() I get:
 
 3: get(dataname)
 2: tslm(xshort ~ trend) at #19
 1: model.cross.validation(l[[MEN]]$series)
 
 Which points to the call to tslm above.  Since I am not supply 'data' to
 the tslm call (in the forecast package),, I am assuming that the code is
 dying here (in tslm):
 
 if (missing(data)) {
 dataname - as.character(formula)[2]
 x - get(dataname)
 data - data.frame(x)
 colnames(data) - dataname
 }
 
 Any ideas what is failing?
 
 Kevin
 
[[alternative HTML version deleted]]
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 
 
 
 -- 
 Joshua Wiley
 Ph.D. Student, Health Psychology
 Programmer Analyst II, ATS Statistical Consulting Group
 University of California, Los Angeles
 https://joshuawiley.com/

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] multivariate modeling codes

2011-11-11 Thread yurirouge
HI, 

I am relatively new to R and would appreciate some help or directions for
this. 
I am trying to model 3 longitudinal outcomes jointly and to identify some
predictors for these 3 joint outcomes (all continuous). I am trying to find
some codes that I may modify to do this but cannot seem to find anything. 

--
View this message in context: 
http://r.789695.n4.nabble.com/multivariate-modeling-codes-tp4032624p4032624.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] Help

2011-11-11 Thread Francesca
Dear Contributors
I would like to perform this operation using a loop, instead of repeating
the same operation many times.
The numbers from 1 to 4 related to different groups that are in the
database and for which I have the same data.


x-c(1,3,7)

datiP1 - datiP[datiP$city ==1,x];

datiP2 - datiP[datiP$city ==2,x];

datiP3 - datiP[datiP$city ==3,x]

datiP4 - datiP[datiP$city ==4,x];
-- 

Thank you for any help you can provide.

Francesca

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] multivariate modeling codes

2011-11-11 Thread B77S
Not sure if this helps, but did you try Google?

http://www.jstatsoft.org/v35/i09/paper

http://cran.r-project.org/web/packages/lcmm/lcmm.pdf

http://www.warwick.ac.uk/statsdept/useR-2011/abstracts/010411-liquetbenoit.pdf



yurirouge wrote:
 
 HI, 
 
 I am relatively new to R and would appreciate some help or directions for
 this. 
 I am trying to model 3 longitudinal outcomes jointly and to identify some
 predictors for these 3 joint outcomes (all continuous). I am trying to
 find some codes that I may modify to do this but cannot seem to find
 anything.
 


--
View this message in context: 
http://r.789695.n4.nabble.com/multivariate-modeling-codes-tp4032624p4032643.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] multivariate modeling codes

2011-11-11 Thread Joshua Wiley
Hi,

It is not at all clear to me what type of model you are trying to fit.
 You could consider ?manova for a multivariate analysis of variance

require(lme4)
?lmer

for longitudinal regression type models

or perhaps

OpenMx: http://openmx.psyc.virginia.edu/

for structural equation modelling in R.

Cheers,

Josh

On Fri, Nov 11, 2011 at 12:25 PM, yurirouge lilysh...@msn.com wrote:
 HI,

 I am relatively new to R and would appreciate some help or directions for
 this.
 I am trying to model 3 longitudinal outcomes jointly and to identify some
 predictors for these 3 joint outcomes (all continuous). I am trying to find
 some codes that I may modify to do this but cannot seem to find anything.

 --
 View this message in context: 
 http://r.789695.n4.nabble.com/multivariate-modeling-codes-tp4032624p4032624.html
 Sent from the R help mailing list archive at Nabble.com.

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.




-- 
Joshua Wiley
Ph.D. Student, Health Psychology
Programmer Analyst II, ATS Statistical Consulting Group
University of California, Los Angeles
https://joshuawiley.com/

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Help

2011-11-11 Thread Joshua Wiley
Hi Francesca,

Try something like this:

x - c(1, 3, 7)
dati - lapply(1:4, function(i) {datiP[datiP$city == i, x]})

dati[[1]] # datiP1
dati[[2]] # datiP2
dati[[3]] # datiP3
dati[[4]] # datiP4

if the *only* groups are 1, 2, 3, 4 (i.e., 1:4 is exhaustive), this
can be simplified:


dati - by(datiP, datiP$city, `[`, x)

For documentation, see:

?lapply
?by
?[ # to see how I use the extraction operator as a function

Hope this helps,

Josh

On Fri, Nov 11, 2011 at 12:24 PM, Francesca
francesca.panco...@gmail.com wrote:
 Dear Contributors
 I would like to perform this operation using a loop, instead of repeating
 the same operation many times.
 The numbers from 1 to 4 related to different groups that are in the
 database and for which I have the same data.


    x-c(1,3,7)

 datiP1 - datiP[datiP$city ==1,x];

 datiP2 - datiP[datiP$city ==2,x];

 datiP3 - datiP[datiP$city ==3,x]

 datiP4 - datiP[datiP$city ==4,x];
 --

 Thank you for any help you can provide.

 Francesca

        [[alternative HTML version deleted]]

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.




-- 
Joshua Wiley
Ph.D. Student, Health Psychology
Programmer Analyst II, ATS Statistical Consulting Group
University of California, Los Angeles
https://joshuawiley.com/

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Predicting x from y

2011-11-11 Thread Ted Harding
On 11-Nov-11 14:51:19, Schreiber, Stefan wrote:
 Dear list members,
 
 I have just a quick question:
 
 I fitted a non-linear model y=a/x+b to describe my data
 (x=temperature and y=damage in %) and it works really nicely
 (see example below). I have 7 different species and 8 individuals
 per species. I measured damage for each individual per species
 at 4 different temperatures (e.g. -5, -10, -20, -40).
 Using the individuals per species, I fitted one model per species.
 Now I'd like to use the fitted model to go back and predict the
 temperature that causes 50% damage (and it's error). Basically,
 it pretty easy by just rearranging the formula to x=a/(y-b).
 But that way I don't get a measure of that temperature's error,
 do I? Can I use the residual standard error that R gave me for
 the non-linear model fit? Or do I have to fit 8 lines (each
 individual) per species, calculate x based on the 8 individuals
 and then take the mean?
 
 Unfortunately, dose.p from the MASS package doesn't work for
 non-linear models. When I take the log(abs(x)) the relationship
 becomes not satisfactory linear either.
 
 Any suggestions are highly appreciated!
 
 Thank you!
 Stefan
 
 EXAMPLE for species #1:
 
 y.damage-c(5.7388985,1.7813519,3.7321461,2.9671031,
 0.3223196,0.3207941,-1.4197658,-5.3472160,
 41.1826677,29.3115243,31.3208841,35.3934115,
 58.5848778,31.1541049,42.2983479,27.0615648,
 64.1037728,54.7003353,66.7317044,65.4725881,
 72.5755056,67.2683495,64.8717942,65.9603073,
 75.0762273,56.7041960,60.0049429,70.0286506,
 73.2801947,72.7015642,75.0944694,81.0361280)
 
 x.temp-c(-5,-5,-5,-5,-5,-5,-5,-5,-10,-10,-10,-10,-10,-10,-10,
 -10,-20,-20,-20,-20,-20,-20,-20,-20,-40,-40,-40,-40,-40,
 -40,-40,-40)
 
 nls(y.damage~a/x.temp+b,start=list(a=400,b=80))
 plot(y.damage~x.temp,xlab='Temperature',ylab='Damage [%]')
 curve(409.61/x+81.84,from=min(x.temp),to=max(x.temp),add=T)

A couple of comments.

First, in general it is not straightforward to estimate
the value of a covariate (here temperature) by inverting
the regression of a response (here damage) on that covariate.
This the inverse regression or calibration problem,
(and it is problematic)! For instance, in linear regression
the estimate obtained by inversion has (theoretically)
no expectation, and has infinite variance. For an outline,
and a few references, see the Wikipedia article:

  http://en.wikipedia.org/wiki/Calibration_(statistics)

Second, I would be inclined to try nls() on a reformulated
version of the problem. Let T50 denote the temperature for
50% damage, and introduce this as a parameter (displacing
your parameter a):

  y = 50*(b + T50)/(b + x)

where T50 = a/50 - b in terms of your original parameters
a and b. With this formula for the non-linear dependence
of damage on temperature, it is no longer necessAry to invert
the regression equation, since the parameter you want is
already there and will be estimated directly.

Hoping this helps,
Ted.


E-Mail: (Ted Harding) ted.hard...@wlandres.net
Fax-to-email: +44 (0)870 094 0861
Date: 11-Nov-11   Time: 21:15:57
-- XFMail --

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] (no subject)

2011-11-11 Thread Simon Kiss
Dear colleagues,
I'm trying to fit a multinomial logistic regression for an ordinal variable.  
I see in the help pages for multinom in nnet that one should scale the 
predictors from 0-1.  Is that really necessary?
Also: can anyone clarify what the difference between alternative-specific and 
individual specific variables are?
Yours, Simon Kiss



l
*
Simon J. Kiss, PhD
Assistant Professor, Wilfrid Laurier University
73 George Street
Brantford, Ontario, Canada
N3T 2C9
Cell: +1 905 746 7606

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] performance of adaptIntegrate vs. integrate

2011-11-11 Thread baptiste auguie
Dear Ravi,

Thank you for your answer.

The integrand I proposed was a dummy example for demonstration
purposes. I experienced a similar slowdown in a real problem, where
knowing in advance the shape of the integrand would not be so easy.

Your advice is sound; I would have to study the underlying code of the
two implementations to know where the difference lies. Delving into
the source code and the algorithms gets quite technical though, so I
was hoping someone already familiar with integrate's internals might
shed some light.

Thanks,

baptiste



On 12 November 2011 03:55, Ravi Varadhan rvarad...@jhmi.edu wrote:
 The integrand is highly peaked.  It is approximately an impulse function
 where much of the mass is concentrated at a very small interval.  Plot the
 function and see for yourself.  This is the likely cause of the problem.



 Other types of integrands where you could experience problems are:
 integrands with singularity at either limit and slowly decaying oscillatory
 integrands.  As to why integrate performs better than adaptIntegrate in this
 situation, I don’t know.  You have to study the two implementations.
  Wynn’s epsilon algorithm is an extrapolation method for improving the
 convergence of a sequence.  This could be an explanation for the better
 performance, but I cannot say for sure.



 Hope this is helpful,

 Ravi

 ---

 Ravi Varadhan, Ph.D.

 Assistant Professor,

 Division of Geriatric Medicine and Gerontology School of Medicine Johns
 Hopkins University



 Ph. (410) 502-2619

 email: rvarad...@jhmi.edu



__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] performance of adaptIntegrate vs. integrate

2011-11-11 Thread baptiste auguie
Dear Hans,

[see inline below]

On 11 November 2011 22:44, Hans W Borchers hwborch...@googlemail.com wrote:
 baptiste auguie baptiste.auguie at googlemail.com writes:


 Dear list,

 [cross-posting from Stack Overflow where this question has remained
 unanswered for two weeks]

 I'd like to perform a numerical integration in one dimension,

 I = int_a^b f(x) dx

 where the integrand f: x in IR - f(x) in IR^p is vector-valued.
 integrate() only allows scalar integrands, thus I would need to call
 it many (p=200 typically) times, which sounds suboptimal. The cubature
 package seems well suited, as illustrated below,

 library(cubature)
 Nmax - 1e3
 tolerance - 1e-4
 integrand - function(x, a=0.01) c(exp(-x^2/a^2), cos(x))
 adaptIntegrate(integrand, -1, 1, tolerance, 2, max=Nmax)$integral
 [1] 0.01772454 1.68294197

 However, adaptIntegrate appears to perform quite poorly when compared
 to integrate. Consider the following example (one-dimensional
 integrand),

 library(cubature)
 integrand - function(x, a=0.01) exp(-x^2/a^2)*cos(x)
 Nmax - 1e3
 tolerance - 1e-4

 # using cubature's adaptIntegrate
 time1 - system.time(replicate(1e3, {
   a - adaptIntegrate(integrand, -1, 1, tolerance, 1, max=Nmax)
 }) )

 # using integrate
 time2 - system.time(replicate(1e3, {
   b - integrate(integrand, -1, 1, rel.tol=tolerance, subdivisions=Nmax)
 }) )

 time1
 user  system elapsed
   2.398   0.004   2.403
 time2
 user  system elapsed
   0.204   0.004   0.208

 a$integral
  [1] 0.0177241
 b$value
  [1] 0.0177241

 a$functionEvaluations
  [1] 345
 b$subdivisions
  [1] 10

 Somehow, adaptIntegrate was using many more function evaluations for a
 similar precision. Both methods apparently use Gauss-Kronrod
 quadrature, though ?integrate adds a Wynn's Epsilon algorithm. Could
 that explain the large timing difference?


 Cubature is astonishingly slow here though it was robust and accurate in
 most cases I used it. You may write to the maintainer for more information.


Indeed, I have been a happy user of cubature too, for
higher-dimensional problems. And I've used other codes from Steve
Johnson before which were all of high quality. I might send an email
to both him and the package developer to inquire about this poor
performance.

 Even a pure R implementation of adaptive Gauss-Kronrod as in pracma::quadgk
 is faster than cubature:

    library(pracma)
    time3 - system.time(replicate(1e3,
      { d - quadgk(integrand, -1, 1) }        # 0.0177240954011303
    ) )
    #  user  system  elapsed
    # 0.899   0.002    0.897

 If your functions are somehow similar and you are willing to dispense with
 the adaptive part, then compute Gaussian quadrature nodes xi and weights wi
 for the fixed interval and compute sum(wi * fj(xi)) for each function fj.
 This avoids recomputing nodes and weights for every call or function.

I've used such a technique before, in C++ code intefaced with Rcpp.
however, I do like the adaptative part, and implementing it seems less
trivial. I don't really want to reinvent the wheel if I can help it
find a faster track.


 Package 'gaussquad' provides such nodes and weights. Or copy the relevant
 calculations from quadgk (the usual (7,15)-rule).

I've used the nodes and weights from statmod::gauss.quad.

Thanks for the tips.

Best regards,

baptiste


 Regards, Hans Werner


 I'm open to suggestions of alternative ways of dealing with
 vector-valued integrands.

 Thanks.

 baptiste



 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Encryption package of R

2011-11-11 Thread Ben Bolker
Jun Ji junj at stanford.edu writes:

 Is there any encryption package of R? For instance, 3DES or AES algorithm
 package or whatever...
 Thanks a lot in advance.
 

  library(sos)
  findFn(encryption)

 finds a partial AES implementation, but not much else.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Win upgrade pb (virus)

2011-11-11 Thread John C Frain
I would suggestto Mario that, if he wishes to use colorspace and keep
AVG,  he should install R 2.13.2.  The colorspace library built with
that version of R does not give a false positive with AVG.  To update
paackages he must remove not only colorspace but also its reverse
dependencies, the reverse dependencies of its dependencies etc.

 I have done this and am running 2.13.2 and 2.14 in parallel.  If you
upgrade R to a new version it is always worthwhile to keep the older
version until you are sure that the new version can run all the
packages that you need.  Apart from problems with anti-virus software
there is a posibility that there may be a delay in  one of your
packages or a dependency becoming available.

I have reported the false positive to AVG and have received an
acknowledgment from them.  The AVG program contains a facility for
reporting false positives. Please use it and put a little more
pressure on them to update their software.

The web site http://virusscan.jotti.org allows you to submit a problem
file to 20 different anti-virus programs.  Only 1 program (AVG) finds
a problem in colorspace.dll.  The other 19 pass the file.

Best Regards

John

On 11 November 2011 10:54, Prof Brian Ripley rip...@stats.ox.ac.uk wrote:
 On Fri, 11 Nov 2011, Mario Valle wrote:

 I just upgraded my Win7 32bits installation to 2.14.0 after deinstalling
 2.12.x
 First thing I moved the win-library from 2.12 to 2.14 and executed a
 update.packages(ask='graphics',checkBuilt=TRUE)
 (Swiss mirror).
 This aborts with the console message:

 Error in if (any(diff)) { : missing value where TRUE/FALSE needed

 And the antivirus (AVG) pops up complaining that colorspace.dll contains
 the virus Win32/Heur
 I cannot ignore the thread because update.packages has already crashed
 when I can reach the antivirus dialogbox to push the ignore button.
 To continue I had to remove the colorspace and ggplot2 packages.

 Is it a known problem? What else can I do (except removing the above
 packages or changing antivirus)?

 Yes (a known problem in AVG), and the standard advice is to change antivirus
 (or set exceptions, which I gather AVG does not allow). And please report
 the false positive to AVG: the more they hear about the more attention they
 will give it.

 The binary distribution has been checked by e.g. Sophos (which both
 winbuilder and ox.ac.uk use).


 Thanks!
                           mario

 --
 Ing. Mario Valle
 Data Analysis and Visualization Group            |
 http://www.cscs.ch/~mvalle
 Swiss National Supercomputing Centre (CSCS)      | Tel:  +41 (91)
 610.82.60
 v. Cantonale Galleria 2, 6928 Manno, Switzerland | Fax:  +41 (91)
 610.82.82


 --
 Brian D. Ripley,                  rip...@stats.ox.ac.uk
 Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
 University of Oxford,             Tel:  +44 1865 272861 (self)
 1 South Parks Road,                     +44 1865 272866 (PA)
 Oxford OX1 3TG, UK                Fax:  +44 1865 272595

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.




-- 
John C Frain
Economics Department
Trinity College Dublin
Dublin 2
Ireland
www.tcd.ie/Economics/staff/frainj/home.html
mailto:fra...@tcd.ie
mailto:fra...@gmail.com

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] Combining Overlapping Data

2011-11-11 Thread kickout
I've scoured the archives but have found no concrete answer to my question.

Problem: Two data sets

1st data set(x) = 20,000 rows 
2nd data set(y) = 5,000 rows

Both have the same column names, the column of interest to me is a variable
called strain.

For example, a strain named Chab1405 appears in x 150 times and in y 25
times...
strain Chab1999 only appears 200 times in x and none in y (so i dont want
that retained).


I want to create a new data frame that has all 175 measurements for
Chab1405 and any other 'strain' that appears in both the two data sets..
but not strains that appear in only one data set...So i want the
intersection of two data sets (maybe?).

I've tried x %in% y, but that only gives TRUE/FALSE


--
View this message in context: 
http://r.789695.n4.nabble.com/Combining-Overlapping-Data-tp4032719p4032719.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] need help

2011-11-11 Thread Supreet kaur
hello all R experts,
 how do I calculate the reliability between the two groups
using the ICCs?

I'll appreciate your reply,

Thanks
Sincerely,

Supreet kaur,
Biomedical research engineer,
Nationwide Childrens Hospital,
Columbus, OH
(614)355-3509

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] list.dir() function

2011-11-11 Thread Mary Kindall
Hi
I have an organism directory that contains two folders galGal3 and hg19 and
many other files.

orgDir = '/home/mary/org'

When I try to use list.dir() function, it gives me the same answer, no
matter what is the value of full.names argument.


 list.dirs(path = indexDir, full.names = FALSE)[1] /home/mary/org

[2] /home/mary/org/galGal3

[3] /home/mary/org/hg19
 list.dirs(path = indexDir, full.names = TRUE)

[1] /home/mary/org

[2] /home/mary/org/galGal3

[3] /home/mary/org/hg19


Also, It prints the directory itself which I don't want to be printed.


Why it is so? Any workaround for this problem?



Thanks


-- 
-
Mary Kindall
Yorktown Heights, NY
USA

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Predicting x from y

2011-11-11 Thread Ted Harding
Follow-up: See at end.

On 11-Nov-11 21:16:02, Ted Harding wrote:
 On 11-Nov-11 14:51:19, Schreiber, Stefan wrote:
 Dear list members,
 
 I have just a quick question:
 
 I fitted a non-linear model y=a/x+b to describe my data
 (x=temperature and y=damage in %) and it works really nicely
 (see example below). I have 7 different species and 8 individuals
 per species. I measured damage for each individual per species
 at 4 different temperatures (e.g. -5, -10, -20, -40).
 Using the individuals per species, I fitted one model per species.
 Now I'd like to use the fitted model to go back and predict the
 temperature that causes 50% damage (and it's error). Basically,
 it pretty easy by just rearranging the formula to x=a/(y-b).
 But that way I don't get a measure of that temperature's error,
 do I? Can I use the residual standard error that R gave me for
 the non-linear model fit? Or do I have to fit 8 lines (each
 individual) per species, calculate x based on the 8 individuals
 and then take the mean?
 
 Unfortunately, dose.p from the MASS package doesn't work for
 non-linear models. When I take the log(abs(x)) the relationship
 becomes not satisfactory linear either.
 
 Any suggestions are highly appreciated!
 
 Thank you!
 Stefan
 
 EXAMPLE for species #1:
 
 y.damage-c(5.7388985,1.7813519,3.7321461,2.9671031,
 0.3223196,0.3207941,-1.4197658,-5.3472160,
 41.1826677,29.3115243,31.3208841,35.3934115,
 58.5848778,31.1541049,42.2983479,27.0615648,
 64.1037728,54.7003353,66.7317044,65.4725881,
 72.5755056,67.2683495,64.8717942,65.9603073,
 75.0762273,56.7041960,60.0049429,70.0286506,
 73.2801947,72.7015642,75.0944694,81.0361280)
 
 x.temp-c(-5,-5,-5,-5,-5,-5,-5,-5,-10,-10,-10,-10,-10,-10,-10,
 -10,-20,-20,-20,-20,-20,-20,-20,-20,-40,-40,-40,-40,-40,
 -40,-40,-40)
 
 nls(y.damage~a/x.temp+b,start=list(a=400,b=80))
 plot(y.damage~x.temp,xlab='Temperature',ylab='Damage [%]')
 curve(409.61/x+81.84,from=min(x.temp),to=max(x.temp),add=T)
 
 A couple of comments.
 
 First, in general it is not straightforward to estimate
 the value of a covariate (here temperature) by inverting
 the regression of a response (here damage) on that covariate.
 This the inverse regression or calibration problem,
 (and it is problematic)! For instance, in linear regression
 the estimate obtained by inversion has (theoretically)
 no expectation, and has infinite variance. For an outline,
 and a few references, see the Wikipedia article:
 
   http://en.wikipedia.org/wiki/Calibration_(statistics)
 
 Second, I would be inclined to try nls() on a reformulated
 version of the problem. Let T50 denote the temperature for
 50% damage, and introduce this as a parameter (displacing
 your parameter a):
 
   y = 50*(b + T50)/(b + x)
 
 where T50 = a/50 - b in terms of your original parameters
 a and b. With this formula for the non-linear dependence
 of damage on temperature, it is no longer necessAry to invert
 the regression equation, since the parameter you want is
 already there and will be estimated directly.
 
 Hoping this helps,
 Ted.

I think I have mis-read your model: I read it as

  y = a/(x+b)

whereas you wrote y/x+b and your model formula in nls()
is y.damage~a/x.temp+b, i.e. y.damage ~ (a/x.temp) + b
which confirms it.

In that case, you may be able to get a satisfactory result
by using a linear regression with

  y.damage = a*z.temp + b

where z.temp = 1/x.temp so the model formula would be

  y.damage ~ z.temp

You then have the straightforward inverse regression
problem (aka calibration problem). The solution to this
takes a bit of explanation, which I do not have the time
for right now. I will write further about it in the morning.

Ted.


E-Mail: (Ted Harding) ted.hard...@wlandres.net
Fax-to-email: +44 (0)870 094 0861
Date: 11-Nov-11   Time: 23:07:35
-- XFMail --

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Combining Overlapping Data

2011-11-11 Thread Sarah Goslee
What about merge() with all=FALSE?

 x - data.frame(a=letters[1:6], b=1:6)
 y - data.frame(a=letters[4:9], b=11:16)
 x
  a b
1 a 1
2 b 2
3 c 3
4 d 4
5 e 5
6 f 6
 y
  a  b
1 d 11
2 e 12
3 f 13
4 g 14
5 h 15
6 i 16
 merge(x, y, by=a, all=FALSE)
  a b.x b.y
1 d   4  11
2 e   5  12
3 f   6  13


If that doesn't work, some sample data would be useful.

Sarah

On Fri, Nov 11, 2011 at 4:07 PM, kickout kyle.ko...@gmail.com wrote:
 I've scoured the archives but have found no concrete answer to my question.

 Problem: Two data sets

 1st data set(x) = 20,000 rows
 2nd data set(y) = 5,000 rows

 Both have the same column names, the column of interest to me is a variable
 called strain.

 For example, a strain named Chab1405 appears in x 150 times and in y 25
 times...
 strain Chab1999 only appears 200 times in x and none in y (so i dont want
 that retained).


 I want to create a new data frame that has all 175 measurements for
 Chab1405 and any other 'strain' that appears in both the two data sets..
 but not strains that appear in only one data set...So i want the
 intersection of two data sets (maybe?).

 I've tried x %in% y, but that only gives TRUE/FALSE


-- 
Sarah Goslee
http://www.functionaldiversity.org

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] need help

2011-11-11 Thread Sarah Goslee
On Fri, Nov 11, 2011 at 4:41 PM, Supreet kaur ksupre...@gmail.com wrote:
 hello all R experts,
                 how do I calculate the reliability between the two groups
 using the ICCs?

Possibly by using GmeanRel() from the multilevel package.

Searching for
reliability group ICC
at http://www.rseek.org
offers several other options, one of which may be more appropriate for
your unspecified application.

Sarah





-- 
Sarah Goslee
http://www.functionaldiversity.org

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Combining Overlapping Data

2011-11-11 Thread jim holtman
Use 'intersect' to get the items common in both dataframes and then use
that to extract the data in common.

On Friday, November 11, 2011, kickout kyle.ko...@gmail.com wrote:
 I've scoured the archives but have found no concrete answer to my
question.

 Problem: Two data sets

 1st data set(x) = 20,000 rows
 2nd data set(y) = 5,000 rows

 Both have the same column names, the column of interest to me is a
variable
 called strain.

 For example, a strain named Chab1405 appears in x 150 times and in y 25
 times...
 strain Chab1999 only appears 200 times in x and none in y (so i dont
want
 that retained).


 I want to create a new data frame that has all 175 measurements for
 Chab1405 and any other 'strain' that appears in both the two data sets..
 but not strains that appear in only one data set...So i want the
 intersection of two data sets (maybe?).

 I've tried x %in% y, but that only gives TRUE/FALSE


 --
 View this message in context:
http://r.789695.n4.nabble.com/Combining-Overlapping-Data-tp4032719p4032719.html
 Sent from the R help mailing list archive at Nabble.com.

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


-- 
Jim Holtman
Data Munger Guru

What is the problem that you are trying to solve?
Tell me what you want to do, not how you want to do it.

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Predicting x from y

2011-11-11 Thread Schreiber, Stefan
Thanks Ted!

I really appreciate your time!

Thanks for the link about the 'problem of calibration', and your suggestion to 
reformulate my model. I had no idea about it before. I certainly learnt 
something today.

I will try your suggestions later today and let you know how it works out.

Stefan


-Original Message-
From: r-help-boun...@r-project.org on behalf of ted.hard...@wlandres.net
Sent: Fri 11/11/2011 4:07 PM
To: r-help@r-project.org
Subject: Re: [R] Predicting x from y
 
Follow-up: See at end.

On 11-Nov-11 21:16:02, Ted Harding wrote:
 On 11-Nov-11 14:51:19, Schreiber, Stefan wrote:
 Dear list members,
 
 I have just a quick question:
 
 I fitted a non-linear model y=a/x+b to describe my data
 (x=temperature and y=damage in %) and it works really nicely
 (see example below). I have 7 different species and 8 individuals
 per species. I measured damage for each individual per species
 at 4 different temperatures (e.g. -5, -10, -20, -40).
 Using the individuals per species, I fitted one model per species.
 Now I'd like to use the fitted model to go back and predict the
 temperature that causes 50% damage (and it's error). Basically,
 it pretty easy by just rearranging the formula to x=a/(y-b).
 But that way I don't get a measure of that temperature's error,
 do I? Can I use the residual standard error that R gave me for
 the non-linear model fit? Or do I have to fit 8 lines (each
 individual) per species, calculate x based on the 8 individuals
 and then take the mean?
 
 Unfortunately, dose.p from the MASS package doesn't work for
 non-linear models. When I take the log(abs(x)) the relationship
 becomes not satisfactory linear either.
 
 Any suggestions are highly appreciated!
 
 Thank you!
 Stefan
 
 EXAMPLE for species #1:
 
 y.damage-c(5.7388985,1.7813519,3.7321461,2.9671031,
 0.3223196,0.3207941,-1.4197658,-5.3472160,
 41.1826677,29.3115243,31.3208841,35.3934115,
 58.5848778,31.1541049,42.2983479,27.0615648,
 64.1037728,54.7003353,66.7317044,65.4725881,
 72.5755056,67.2683495,64.8717942,65.9603073,
 75.0762273,56.7041960,60.0049429,70.0286506,
 73.2801947,72.7015642,75.0944694,81.0361280)
 
 x.temp-c(-5,-5,-5,-5,-5,-5,-5,-5,-10,-10,-10,-10,-10,-10,-10,
 -10,-20,-20,-20,-20,-20,-20,-20,-20,-40,-40,-40,-40,-40,
 -40,-40,-40)
 
 nls(y.damage~a/x.temp+b,start=list(a=400,b=80))
 plot(y.damage~x.temp,xlab='Temperature',ylab='Damage [%]')
 curve(409.61/x+81.84,from=min(x.temp),to=max(x.temp),add=T)
 
 A couple of comments.
 
 First, in general it is not straightforward to estimate
 the value of a covariate (here temperature) by inverting
 the regression of a response (here damage) on that covariate.
 This the inverse regression or calibration problem,
 (and it is problematic)! For instance, in linear regression
 the estimate obtained by inversion has (theoretically)
 no expectation, and has infinite variance. For an outline,
 and a few references, see the Wikipedia article:
 
   http://en.wikipedia.org/wiki/Calibration_(statistics)
 
 Second, I would be inclined to try nls() on a reformulated
 version of the problem. Let T50 denote the temperature for
 50% damage, and introduce this as a parameter (displacing
 your parameter a):
 
   y = 50*(b + T50)/(b + x)
 
 where T50 = a/50 - b in terms of your original parameters
 a and b. With this formula for the non-linear dependence
 of damage on temperature, it is no longer necessAry to invert
 the regression equation, since the parameter you want is
 already there and will be estimated directly.
 
 Hoping this helps,
 Ted.

I think I have mis-read your model: I read it as

  y = a/(x+b)

whereas you wrote y/x+b and your model formula in nls()
is y.damage~a/x.temp+b, i.e. y.damage ~ (a/x.temp) + b
which confirms it.

In that case, you may be able to get a satisfactory result
by using a linear regression with

  y.damage = a*z.temp + b

where z.temp = 1/x.temp so the model formula would be

  y.damage ~ z.temp

You then have the straightforward inverse regression
problem (aka calibration problem). The solution to this
takes a bit of explanation, which I do not have the time
for right now. I will write further about it in the morning.

Ted.


E-Mail: (Ted Harding) ted.hard...@wlandres.net
Fax-to-email: +44 (0)870 094 0861
Date: 11-Nov-11   Time: 23:07:35
-- XFMail --

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide 

Re: [R] multivariate modeling codes

2011-11-11 Thread Dennis Murphy
Hi:

R-Bloggers picked up on this web site that contains several
interesting posts re usage of R, including this one on using a
Bayesian approach to multivariate mixed effects models using the
MCMCglmm package:
http://www.quantumforest.com/2011/11/coming-out-of-the-bayesian-closet-multivariate-version/

With multiple responses in a longitudinal study, this might be of
interest to you. A similar question re multivariate response was asked
a few days ago on R-sig-mixed-models, where Kevin Wright posted the
above link. I'd suggest that future questions on this subject be moved
to that group.

HTH,
Dennis

On Fri, Nov 11, 2011 at 12:25 PM, yurirouge lilysh...@msn.com wrote:
 HI,

 I am relatively new to R and would appreciate some help or directions for
 this.
 I am trying to model 3 longitudinal outcomes jointly and to identify some
 predictors for these 3 joint outcomes (all continuous). I am trying to find
 some codes that I may modify to do this but cannot seem to find anything.

 --
 View this message in context: 
 http://r.789695.n4.nabble.com/multivariate-modeling-codes-tp4032624p4032624.html
 Sent from the R help mailing list archive at Nabble.com.

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Combining Overlapping Data

2011-11-11 Thread Dennis Murphy
Hi:

This doesn't sort the data by strain level, but I think it does what
you're after. It helps if strain is either a factor or character
vector in each data frame.

h - function(x, y) {
   tbx - table(x$strain)
   tby - table(y$strain)
  # Select the strains who have more than one member
  # in each data frame
   mgrps - intersect(names(tbx[tbx  0]),
  names(tby[tby  0]))
  # concatenate the data with common strains
   rbind(subset(x, gp %in% mgrps),
 subset(y, gp %in% mgrps))
   }

# Result:
dc - h(x, y)

HTH,
Dennis

On Fri, Nov 11, 2011 at 1:07 PM, kickout kyle.ko...@gmail.com wrote:
 I've scoured the archives but have found no concrete answer to my question.

 Problem: Two data sets

 1st data set(x) = 20,000 rows
 2nd data set(y) = 5,000 rows

 Both have the same column names, the column of interest to me is a variable
 called strain.

 For example, a strain named Chab1405 appears in x 150 times and in y 25
 times...
 strain Chab1999 only appears 200 times in x and none in y (so i dont want
 that retained).


 I want to create a new data frame that has all 175 measurements for
 Chab1405 and any other 'strain' that appears in both the two data sets..
 but not strains that appear in only one data set...So i want the
 intersection of two data sets (maybe?).

 I've tried x %in% y, but that only gives TRUE/FALSE


 --
 View this message in context: 
 http://r.789695.n4.nabble.com/Combining-Overlapping-Data-tp4032719p4032719.html
 Sent from the R help mailing list archive at Nabble.com.

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] Time series analysis with random effects

2011-11-11 Thread Schrabauke
Hi folks,
i have some problems with my evaluation. We have collect tons of data from
23 testpersons for our new road study.
I have now a time series for each person and all the logs when he
accelerates or hits the break trying to solve five different tasks.

The dataset lools like:
#Unixtime   time   accelerate  break UserID
task 
#10372 1312358742 10:05:41.600 3.054   0.000  2 lane
#10373 1312358742 10:05:41.700 3.056   1.000  2
emergency

The break variable is binary and i want to test if there is a significant
difference between the tasks using the break variable.

I have not really a clue what i have to look for. My idea was to use a 2
Sample t-test for the accelerate variable but the conditions are not
fulfilled for this test.

I am not really experienced with time series analysis so thanks a lot for
your advice!
Bernd





   

--
View this message in context: 
http://r.789695.n4.nabble.com/Time-series-analysis-with-random-effects-tp4033263p4033263.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] predictive apriori

2011-11-11 Thread Flavio Barros
Dear list members,

I know that there is the arules package with the implementation of the
apriori algorithm. However i want to use the predictive apriori instead.
These algorithm can mine as rules as i want and there is an implementation
on weka.

There is some implementation on R?

-- 
Att,

Flávio Barros

[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Encryption package of R

2011-11-11 Thread Marc Schwartz

On Nov 11, 2011, at 4:04 PM, Ben Bolker wrote:

 Jun Ji junj at stanford.edu writes:
 
 Is there any encryption package of R? For instance, 3DES or AES algorithm
 package or whatever...
 Thanks a lot in advance.
 
 
  library(sos)
  findFn(encryption)
 
 finds a partial AES implementation, but not much else.

Hi,

When it comes to encryption, you are best to use a library that has a good 
history behind it.

I would suggest looking at http://www.gnupg.org/ and perhaps specifically, 
libgcrypt (http://directory.fsf.org/wiki/Libgcrypt).

HTH,

Marc Schwartz

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] Odds ratios from lrm plot

2011-11-11 Thread RD235
The code

library(Design)
f - lrm(y~x1+x2+x1*x2, data=data)
plot(f)

produces a plot of log odds vs x2 with 0.95 confidence intervals. How do I
get a plot of odds ratios vs x2 instead?

Thanks


--
View this message in context: 
http://r.789695.n4.nabble.com/Odds-ratios-from-lrm-plot-tp4033340p4033340.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Odds ratios from lrm plot

2011-11-11 Thread Frank Harrell
The Design package is obsolete.  Use its replacement rms - see
http://biostat.mc.vanderbilt.edu/Rrms.  Use something like the following to
reproduce what you already have:

To plot odds ratios against a specific value of x2 and at a specific x1,
type ?contrast.rms to see examples.
Frank

RD235 wrote:
 
 The code
 
 library(Design)
 f - lrm(y~x1+x2+x1*x2, data=data)
 plot(f)
 
 produces a plot of log odds vs x2 with 0.95 confidence intervals. How do I
 get a plot of odds ratios vs x2 instead?
 
 Thanks
 


-
Frank Harrell
Department of Biostatistics, Vanderbilt University
--
View this message in context: 
http://r.789695.n4.nabble.com/Odds-ratios-from-lrm-plot-tp4033340p4033414.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Odds ratios from lrm plot

2011-11-11 Thread David Winsemius


On Nov 11, 2011, at 9:08 PM, RD235 wrote:


The code

library(Design)
f - lrm(y~x1+x2+x1*x2, data=data)
plot(f)

produces a plot of log odds vs x2 with 0.95 confidence intervals.  
How do I

get a plot of odds ratios vs x2 instead?


You would construct a dataset that had selected combinations of x1 and  
x2 and then plot the response.


It's been a couple of years since I used Design. Harrell migrated to  
lattice plots when he upgraded to 'rms' and changed the Predict/plot  
interface a bit. In the new version it might  something like  
plot(predict(f, x1=seq( ...), x2=c(1,2), type= fitted), deopending  
on teh data arrangement.  Why don't you look more closely at the help  
page for `lrm` and work through the examples. At least in the current  
'rms' version there is more than one worked example, and I'm pretty  
sure there were one in Design before that.



--


David Winsemius, MD
West Hartford, CT

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] fitting a distribtution to event data

2011-11-11 Thread manjax
Hi all

I am trying to fit a distribution (i.e. gamma) to some data and I understand
how to use the fitdistr  from the MASS package. However, what should you do
when the data has a probability associated with it that are not all equal.

e.g.
Pr(Occurrence)  Size
0.00460
0.02 30
0.00670
0.06 40


regards

--
View this message in context: 
http://r.789695.n4.nabble.com/fitting-a-distribtution-to-event-data-tp4033760p4033760.html
Sent from the R help mailing list archive at Nabble.com.
[[alternative HTML version deleted]]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Odp: Error in matrix, not ordered vectors or numerical value, and SIAR.

2011-11-11 Thread Alex DJ
Hello Petr, 

The demo's don't run either, with the same errors. 

Thanks for your help. 

Best wishes. 

--
View this message in context: 
http://r.789695.n4.nabble.com/Error-in-matrix-not-ordered-vectors-or-numerical-value-and-SIAR-tp4024578p4033682.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] lm: how are polynomial functions interpreted?

2011-11-11 Thread geeknick
A 
http://polynomial-trimonial-binomial.blogspot.com/2011/10/how-to-deal-with-polynomials.html
polynomial function  is evaluated by foiling out the equation so you can
solve for 'x'.

--
View this message in context: 
http://r.789695.n4.nabble.com/lm-how-are-polynomial-functions-interpreted-tp875663p4033696.html
Sent from the R help mailing list archive at Nabble.com.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.