Re: [R] Need help to estimate the Coef matrices in mAr

2006-09-04 Thread Spencer Graves
  Have you tried 'RSiteSearch(multivariate autoregression, 
functions)'?  This produced 14 hits for me just now, the first of 
which mentions a package 'MSBVAR'.  Have you looked at that? 

  If that failed, I don't think it would be too hard to modify 
'mAr.est' to do what you want.  If it were my problem, I might a local 
copy of the function, then add an argument accepting a 2 or 
3-dimensional array with numbers for AR coefficients to be fixed and NAs 
for the coefficients.  Then I'd use 'debug' to walk through the function 
line by line until I figured out how to modify the function to do what I 
wanted.  I haven't checked all the details, so I don't know for sure if 
this would work, but the function contains a line 'R = qr.R(qr((rbind(K, 
diag(scale, complete = TRUE)' which I would start by decomposing, 
possibly starting as follows: 

  Z - rbind(K, diag(scale)

I'd figure out how the different columns of Z relate to my problem, then 
modify it appropriately to get what I wanted. 

  Another alternative would be to program it from scratch using 
something like 'optim' to minimize the sum of squares of residuals over 
the free parameters in my AR matrices.   I'm confident I could make this 
work, even if the I somehow could not get it with either of the other two. 

  There may be something else  better, e.g., a Kalman filter 
representation, but I can't think how to do that off the top if my head. 

  Hope this helps. 
  Spencer Graves

Arun Kumar Saha wrote:
 Dear R users,

 I am using mAr package to fit a Vector autoregressive model to my data. But
 here I want to put some predetermined values for some elements in
 coefficient matrix that mAr.est going to estimate. For example if p=3 then I
 want to put A3[1,3] = 0 and keep rest of the elements of coefficient
 matrices to be determined by mAr.est.

 Can anyone please tell me how can I do that?

 Sincerely yours,
 Arun

   [[alternative HTML version deleted]]

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


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


[R] Subsetting vectors based on condition

2006-09-04 Thread Mahesh Krishnan
Hello,

I have a question regarding subsetting of vectors. Here's an example of
what I'm trying to do:

vect.1 - c(76,195, 290, 380)

vect.2 -  c(63,  95, 133, 170, 215, 253, 285, 299, 325, 375)

I would like to subset vect.2 so that it has the same length as vect.1,
and its numbers are the first corresponging higher value compared to
vect.1.

The output should be:

 final.output = (95, 215, 299, NA)

What is the fastest/most eficient way to accompllish this in R?

Thanks for the help in advance,

Mahesh Krishnan

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


[R] Command line cut-off

2006-09-04 Thread Björn Thalheim
Hi,

I noticed that, when working with R on a command line, I cannot enter
anything into the command line which is longer than 1023 characters.
If I input sth longer, it'l be just cut off at this lenght, and a +
will be presented to me on the command line.

It's not as bad since I can copy'n'paste commands longer than this with
newline characters in between onto the prompt, but still kind of annoying.

Does anybody know how I can make my command line accept lines longer
than 1023 characters?

Regards,

Björn


-- 
Q:  How many lawyers does it take to change a light bulb?
A:  One.  Only it's his light bulb when he's done.

-- 
Important! Please recognize my new GPG Public Key!
 Björn Thalheim
gpg fingerprint: 2F22 AAEB 1818 1548 EC78  1AE8 9D2E FCB4 0980 28CC
   download key: wget http://www.ifsr.de/~bjoern/gpg/public_key.asc
   See also: http://www.ifsr.de/~bjoern/gpg/key.html



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


[R] [R-pkgs] New forecasting bundle of packages

2006-09-04 Thread Rob Hyndman
v1.0 of the forecasting bundle of packages is now on CRAN and will 
propagate to mirrors shortly.

The forecasting bundle of R packages provides new forecasting methods, 
and graphical tools for displaying and analysing forecasts. It comprises 
the following packages:

 * forecast: Functions and methods for forecasting.

 * fma: All data sets from Makridakis, Wheelwright and Hyndman 
(1998) Forecasting: methods and applications, Wiley  Sons: New 
York.

 * Mcomp: All data from the M1 and M3 forecast competitions.

Key features:

* a forecast method and class which can be applied to Arima, StructTS, 
HoltWinters and other time series models. This is preferred to predict() 
as it provides output in a consistent format (the forecast class) that 
can be used by other functions.

* automatic univariate time series forecasting based on exponential 
smoothing state space models. This is much more general and flexible 
than HoltWinters().

* automatic ARIMA forecasting based on minimizing the AIC or BIC.

* several new forecasting methods and time series graphics.

Some features of the forecast package were the subject of my talk at 
UseR! in Vienna in June. Slides of the talk are at 
http://www.robhyndman.info/talks/Hyndman_UseR.pdf

Anyone who has been using earlier versions of the packages from my web 
pages should check out the list of changes at 
http://www.robhyndman.info/Rlibrary/forecast/

__
Professor Rob J Hyndman
Department of Econometrics  Business Statistics,
Monash University, VIC 3800, Australia
http://www.robhyndman.info/

___
R-packages mailing list
R-packages@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-packages

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


Re: [R] generating loglogistic distribution in R

2006-09-04 Thread nmi13
Hi Dear,

Can someone please inform me on genreating random variables of Loglogistic, 
and PERT beta distributions?

Thanks for your time and help, in advance.

Regards,
Murthy.

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


Re: [R] Command line cut-off

2006-09-04 Thread Prof Brian Ripley
On Mon, 4 Sep 2006, Björn Thalheim wrote:

 Hi,
 
 I noticed that, when working with R on a command line, I cannot enter
 anything into the command line which is longer than 1023 characters.

Actually, 1000 or 1022 or 1023 bytes, and it depends on your 'command 
line' (and you have not even told us your OS).

You will find details in the R-devel list archives, e.g.

https://stat.ethz.ch/pipermail/r-devel/2006-August/038985.html

 If I input sth longer, it'l be just cut off at this lenght, and a +
 will be presented to me on the command line.
 
 It's not as bad since I can copy'n'paste commands longer than this with
 newline characters in between onto the prompt, but still kind of annoying.

How do you think thousands users of R for a decade have managed?
What are you doing that needs 'commands' (R has function calls) longer 
than 1000 or so characters that need to be on one line?

 Does anybody know how I can make my command line accept lines longer
 than 1023 characters?

See the R-devel version of R for a full description: on some consoles it 
works there.

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


[R] xlsReadWrite 1.0

2006-09-04 Thread Hans-Peter
The first version of xlsReadWrite has been uploaded to CRAN.


-- WHAT IS IT? / R COMMANDS

The package allows you to read and write Excel files natively. The supported
Format is BIFF 8, i.e. Excel from version 97 up to 2003.

• read.xls( file, colNames = TRUE, sheet = 1, type = data.frame, from = 1
)
• write.xls( x, file, colNames = TRUE, sheet = 1, from = 1 )


-- PLATFORM / LICENSE

Currently the package is only available on Windows (see also below).

The license is GPLv2 with an exception to allow the use of a third party
library (flexcel). With this exception the package as a whole becomes a
non-free package (according to the R license non-free packages are not
encouraged but at least possible/tolerated). See also README.

-- TECHNICAL BACKGROUND / PRO VERSION:

The core of this package is written in Delphi (PASCAL) and uses the R
headerfiles for Delphi which I have written earlier (see website below). For
the actual read/write operation in excel a proprietary library (Flexcel) is
used as I wouldn't be able to write this by myself and because for pascal
there are not open source libraries available. (AFAIK. - I am aware of
Apache POI (Java), a library for OpenOffice (C, too much work to translate
to pascal) and also some older one which don't support BIFF8).

I intend to make available (in 2 - 3 weeks) a pro version which has some
more functionality and formal support. As an independant developer I am
certainly happy if people will check this version out or consider a
contribution, but it is not my idea to play silly shareware games, the
open source version works well and I am commited to support it as far as I
can.

-- OTHER PLATFORMS (Linux/Mac)

Technically xlsReadWrite could be crossplatform. Flexcel already has a Kylix
version, which means that it could be made to run on Linux almost
immediately. The IMO better way would be to port the code to *FreePascal*:
this fits better in the open source world of R and this also would make it
possible to supply a Mac version. - While I started to port it, I didn't
succeed in the 2 days I have given me. - If there are any pascal-lovers
out there I would be very interested to get in touch, it would be more fun
to do it together than to try all by myself.


-- MORE INFOS / CONTACT:

At http://treetron.googlepages.com/ you will find an additional file with
testdata/-scripts and more infos (contact, bugs (currently nothing),
suggestions and todos).

I hope you enjoy using this package (I certainly do...) and am looking
forward to any feedback you might have.

Best regards,
Hans-Peter

[[alternative HTML version deleted]]

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


[R] merge files after cor.test

2006-09-04 Thread jia ding
Dear All,

Suppose I have 2 files:
# first one : testid.csv
A
B
C
D
E

 (id- read.table (testid.csv,col.name=c(id)))
  id
1  A
2  B
3  C
4  D
5  E

# second file is the result file I calculate from cor.text, which shows the
correlation coefficient.
 cor.value.t
  1  2   3   4  5
1 1.000  0.2156213  0.31000492  0.22282154  0.1822277
2 0.2156213  1.000 -0.31183893  0.42681488  0.3421535
3 0.3100049 -0.3118389  1. -0.02801885 -0.1307732
4 0.2228215  0.4268149 -0.02801885  1.  0.1454049
5 0.1822277  0.3421535 -0.13077317  0.14540493  1.000

# I tried to merge these two files together.
# what I expected is like this:
A 1.000  0.2156213  0.31000492  0.22282154  0.1822277
B 0.2156213  1.000 -0.31183893  0.42681488  0.3421535
C 0.3100049 -0.3118389  1. -0.02801885 -0.1307732
D 0.2228215  0.4268149 -0.02801885  1.  0.1454049
E 0.1822277  0.3421535 -0.13077317  0.14540493  1.000


# but after I used (output-merge(id,cor.value.t)), which shows 25 lines
(below):
   id 1  2   3   4  5
1   A 1.000  0.2156213  0.31000492  0.22282154  0.1822277
2   B 1.000  0.2156213  0.31000492  0.22282154  0.1822277
3   C 1.000  0.2156213  0.31000492  0.22282154  0.1822277
4   D 1.000  0.2156213  0.31000492  0.22282154  0.1822277
5   E 1.000  0.2156213  0.31000492  0.22282154  0.1822277
6   A 0.2156213  1.000 -0.31183893  0.42681488  0.3421535
7   B 0.2156213  1.000 -0.31183893  0.42681488  0.3421535
8   C 0.2156213  1.000 -0.31183893  0.42681488  0.3421535
9   D 0.2156213  1.000 -0.31183893  0.42681488  0.3421535
10  E 0.2156213  1.000 -0.31183893  0.42681488  0.3421535
11  A 0.3100049 -0.3118389  1. -0.02801885 -0.1307732
12  B 0.3100049 -0.3118389  1. -0.02801885 -0.1307732
13  C 0.3100049 -0.3118389  1. -0.02801885 -0.1307732
14  D 0.3100049 -0.3118389  1. -0.02801885 -0.1307732
15  E 0.3100049 -0.3118389  1. -0.02801885 -0.1307732
16  A 0.2228215  0.4268149 -0.02801885  1.  0.1454049
17  B 0.2228215  0.4268149 -0.02801885  1.  0.1454049
18  C 0.2228215  0.4268149 -0.02801885  1.  0.1454049
19  D 0.2228215  0.4268149 -0.02801885  1.  0.1454049
20  E 0.2228215  0.4268149 -0.02801885  1.  0.1454049
21  A 0.1822277  0.3421535 -0.13077317  0.14540493  1.000
22  B 0.1822277  0.3421535 -0.13077317  0.14540493  1.000
23  C 0.1822277  0.3421535 -0.13077317  0.14540493  1.000
24  D 0.1822277  0.3421535 -0.13077317  0.14540493  1.000
25  E 0.1822277  0.3421535 -0.13077317  0.14540493  1.000


Is anybody know why it outputs 5 times?

Thanks a lot!

Nina

[[alternative HTML version deleted]]

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


Re: [R] my error with augPred

2006-09-04 Thread Petr Pikal
Hallo

thank you for your response. I am not sure but maybe fixed effects 
cannot be set to be influenced by a factor to be able to use augPred.

lob-Loblolly[Loblolly$Seed!=321,]
set.seed(1)
lob-data.frame(lob, x1=sample(letters[1:3], replace=T)) # add a 
#factor
lob-groupedData(height~age|Seed, data=lob)
fm1 - nlme(height ~ SSasymp(age, Asym, R0, lrc),
data = lob,
fixed = Asym + R0 + lrc ~ 1,
random = Asym ~ 1,
start = c(Asym = 103, R0 = -8.5, lrc = -3.3))

fm2-update(fm1, fixed=list(Asym~x1, R0+lrc~1), start=c(103,0,-8.5,-
3))
 ^^^
and

plot(augPred(fm2))

Throws an error.
So it is not possible to use augPred with such constructions.

Best regards.
Petr Pikal

On 2 Sep 2006 at 17:58, Spencer Graves wrote:

Date sent:  Sat, 02 Sep 2006 17:58:05 -0700
From:   Spencer Graves [EMAIL PROTECTED]
To: Petr Pikal [EMAIL PROTECTED]
Copies to:  r-help@stat.math.ethz.ch
Subject:Re: [R] my error with augPred

 comments in line 
 
 Petr Pikal wrote:
  Dear all
 
  I try to refine my nlme models and with partial success. The model
  is refined and fitted (using Pinheiro/Bates book as a tutorial) but
  when I try to plot
 
  plot(augPred(fit4))
 
  I obtain
  Error in predict.nlme(object, value[1:(nrow(value)/nL), , drop =
  FALSE],  : 
  Levels (0,3.5],(3.5,5],(5,7],(7,Inf] not allowed for 
  vykon.fac

 
  Is it due to the fact that I have unbalanced design with not all
  levels of vykon.fac present in all levels of other explanatory
  factor variable?

 I don't know, but I'm skeptical. 
  I try to repeat 8.19 fig which is OK until I try:
 
  fit4 - update(fit2, fixed = list(A+B~1,xmid~vykon.fac, scal~1), 
  start = c(57, 100, 700, rep(0,3), 13))
 
  I know I should provide an example but maybe somebody will be clever
  enough to point me to an explanation without it.

 I'm not. 
 
 To answer these questions without an example from you, I'd have to
 make up my own example and try to see if I could replicate the error
 messages you report, and I'm not sufficiently concerned about this
 right now to do that. 
 
 Have you tried taking an example from the book and deleting certain
 rows from the data to see if you can force it to reproduce your error?
 
 
 Alternatively, have you tried using 'debug' to trace through the code
 line by line until you learn enough of what it's doing to answer your
 question? 
 
 Spencer Graves
  nlme version 3.1-75
  SSfpl model
  R 2.4.0dev (but is the same in 2.3.1), W2000.
 
  Thank you
  Best regards.
 
  Petr PikalPetr Pikal
  [EMAIL PROTECTED]
 
  __
  R-help@stat.math.ethz.ch mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide
  http://www.R-project.org/posting-guide.html and provide commented,
  minimal, self-contained, reproducible code.
 
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html and provide commented,
 minimal, self-contained, reproducible code.

Petr Pikal
[EMAIL PROTECTED]

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


[R] abline and plot(augPred) help

2006-09-04 Thread Petr Pikal
Dear all

as I did not get any response on my post about abline and 
plot(augPred)) I try again. I hope I do not break some posting guide 
rules. I would try to contact package maintainer directly but there 
is stated to be R-core people, so I feel R-help list shall be OK.

I need to draw straight lines through augPred plotted panels 
(vertical or horizontal) at specified point. I know I shall probably 
use panel.abline but I am missing correct syntax. Below you can see 
my attempts together with results. I hope somebody can point me to 
right direction.

I am probably somewhere close but I have no clue, which parameter I
shall modify to get measured points, fitted lines and vertical lines
in panels together.

Please help

Thank you
Best regards.
Petr Pikal

fm1 - lme(Orthodont)

# standard plot
plot(augPred(fm1, level = 0:1, length.out = 2))

#plot with vertical but without points and fitted lines
plot(augPred(fm1, level = 0:1, length.out = 2),
panel=function(v,...) {
panel.abline(v=10)}
)

# plot with vertical but without fitted lines
plot(augPred(fm1, level = 0:1, length.out=2),
panel=function(x,y,...) {
panel.xyplot(x,y,...)
panel.abline(v=10)}
)

# plot with vertical and with all points (fitted lines are drawn as 
points)
plot(augPred(fm1, level = 0:1),
panel=function(x,y,...) {
panel.xyplot(x,y,...)
panel.abline(v=10)}
)

Petr Pikal
[EMAIL PROTECTED]

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


Re: [R] Subsetting vectors based on condition

2006-09-04 Thread jim holtman
 vect.1
[1]  76 195 290 380
 vect.2
 [1]  63  95 133 170 215 253 285 299 325 375
 unlist(lapply(vect.1, function(x)vect.2[which(vect.2  x)[1]]))
[1]  95 215 299  NA



On 9/3/06, Mahesh Krishnan [EMAIL PROTECTED] wrote:
 Hello,

 I have a question regarding subsetting of vectors. Here's an example of
 what I'm trying to do:

 vect.1 - c(76,195, 290, 380)

 vect.2 -  c(63,  95, 133, 170, 215, 253, 285, 299, 325, 375)

 I would like to subset vect.2 so that it has the same length as vect.1,
 and its numbers are the first corresponging higher value compared to
 vect.1.

 The output should be:

  final.output = (95, 215, 299, NA)

 What is the fastest/most eficient way to accompllish this in R?

 Thanks for the help in advance,

 Mahesh Krishnan

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



-- 
Jim Holtman
Cincinnati, OH
+1 513 646 9390

What is the problem you are trying to solve?

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


Re: [R] Fitting Pareto distribution to some data

2006-09-04 Thread Paul Smith
On 9/3/06, Prof Brian Ripley [EMAIL PROTECTED] wrote:
  I am trying to fit Pareto distribution to some data. MASS package does
  not support Pareto distribution. Is there some alternative way?

 Actually fitdistr{MASS} does if you supply the pdf for a Pareto.
 That is not in base R, but easy to write for yourself.

 It seems that Pareto and generalized Pareto is in several packages,
 including POT SoPhy VaR evd evir fExtremes lmomco.

Thanks, Prof Brian Ripley.

Paul

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


Re: [R] Question on Chi-square of null model in sem package

2006-09-04 Thread John Fox
Dear Wei-Wei,

As I explained to you in private email yesterday (perhaps you didn't receive
my reply?), the problem that you point out is due to a bug in the sem
function that I fixed some time ago and then inadvertently reintroduced.
Yesterday, I sent a corrected version of the sem package (0.9-5) to CRAN;
the source package is there now and I'm sure that the compiled Windows
package will appear in due course.

Thank you once more for bringing the problem to my attention.

John


John Fox
Department of Sociology
McMaster University
Hamilton, Ontario
Canada L8S 4M4
905-525-9140x23604
http://socserv.mcmaster.ca/jfox 
 

 -Original Message-
 From: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] On Behalf Of Guo Wei-Wei
 Sent: Monday, September 04, 2006 12:34 AM
 To: r-help@stat.math.ethz.ch
 Subject: [R] Question on Chi-square of null model in sem package
 
 Dear all,
 
 I met a problem while doing SEM by sem package. I got a 
 negative chi-square of null model. Because the theoretical 
 value of chi-square cannot be negative, I checked the source 
 code of sem.R in sem package and I found the Chi-square of 
 null model was computed by the following
 expression:
 
 result$chisqNull - (N - 1) * (sum(diag(S %*% diag(1/diag(S +
 log(prod(diag(S
 
 I think the reason for negative Chi-square is the too small value of
 prod(diag(S)) of my data. I'm working on a data.frame named 
 emc.data from a sample of a 16-item questioinnaire. The 
 variance of items are
 
  diag(cov(emc.data))
  EMC1  EMC2  EMC3  EMC4  EMC5  EMC6   
EMC7  EMC8
 0.364 0.2350041 0.2488009 0.2901653 0.3195399 0.3107343 
 0.3436622 0.2345912
  EMC9 EMC10 EMC11 EMC12 EMC13 EMC14   
   EMC15 EMC16
 0.2621680 0.3230400 0.4039245 0.3803105 0.2773370 0.4348342 
 0.2757216 0.3405252
 
 The fit indices of RMSEA and GFI are good, so I think the 
 problem might be solve by another way for computing the 
 Chi-square of null model. I'm not well trained in maths, so I 
 come for help. Any advise is appreciated.
 
 Best wishes,
 Wei-Wei
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide 
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

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


[R] library(plgem)

2006-09-04 Thread Amir Safari



 
  Dear Users,
   
  library(plgem) doesn't exist directly in the list of available packages of R. 
Where could it be found?
  Thanks so much for help.
  Amir



-

[[alternative HTML version deleted]]

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


[R] how to fit gauss beam?

2006-09-04 Thread Thomas Hunger
Hello,

I am having a hard time fitting a gauss beam using R. In 
gnutplot I did something like

$ w(z) = w0 * sqrt(1+(z/z0)**2)
$ fit w(z) 'before_eom.txt' using 1:2 via w0, z0

to obtain w0 and z0. Now I want to do the same in R. I tried 
a linear model like this (r = radius, z = distance):

beam - function(z) {
  sum(sqrt(1 + z**2))
}

lm(r ~ I(beam(z)), data = before_eom)

Which gives nonsensical answers ...

Then I tried a nonlinear model:

d - read.table (before_eom.tab, header=T)
z - d$d
r - d$minor * 1e-6

beam- function(p) {
  M - 1.1
  sum((r-M*p[1]*sqrt(1 + (z/p[2])**2))^2)
}

out - nlm(beam, p=c(400, 0.1), hessian=TRUE)
out$estimate

This is very sensitive to the starting values and gives 
nonsensical answers as well...

Is there a simple way to translate the gnuplot fit into R?

Thanks,
Thomas

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


Re: [R] library(plgem)

2006-09-04 Thread john seers \(IFR\)


Perhaps here?


http://www.bioconductor.org/packages/bioc/1.6/src/contrib/html/plgem.htm
l


 

-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Amir Safari
Sent: 04 September 2006 14:44
To: R-help@stat.math.ethz.ch
Subject: [R] library(plgem)





 
  Dear Users,
   
  library(plgem) doesn't exist directly in the list of available
packages of R. Where could it be found?
  Thanks so much for help.
  Amir



-

[[alternative HTML version deleted]]

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

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


Re: [R] newbie question about index

2006-09-04 Thread sun
Thanks for all these replies, all work perfectly.

Sun
--- Petr Pikal [EMAIL PROTECTED]写道:

 Hallo
 
 probably there are other options but
 
 outer(1:3,a, ==)*1
 
 can do what you want.
 
 HTH
 Petr
 
 
 
 On 31 Aug 2006 at 22:41, z s wrote:
 
 Date sent:Thu, 31 Aug 2006 22:41:27 +0800
 (CST)
 From: z s [EMAIL PROTECTED]
 To:   r-help@stat.math.ethz.ch
 Subject:  [R] newbie question about index
 
  Hi,
 
I am trying to convert a variable a =
 sample(1:3,100,rep = T)
represents choices into a 3X100 dummy varible b
 with corresponding
element set to 1 otherwise 0.
  eg.
 
  a: 1 3 2 1 2 3 1 1
 
  b: 1 0 0 1 0 0 1 1..
  0 0 1 0 1 0 0 0...
  0 1 0 0 0 1 0 0...
 
   Is there something like b[a] =1 existing? I could
 not figure this out
   myself.
 
 
  -
   Mp3�杩袼�-新歌热歌高速下
   [[alternative HTML version deleted]]
 
 
 
 Petr Pikal
 [EMAIL PROTECTED]
 
 




___ 
抢注雅虎免费邮箱-3.5G容量,20M附件!

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


Re: [R] generating loglogistic distribution in R

2006-09-04 Thread Ben Bolker
nmi13 nmi13 at ext.canterbury.ac.nz writes:

 
 Hi Dear

   ???

 
 Can someone please inform me on genreating random variables of Loglogistic, 
 and PERT beta distributions?
 

loglogistic:

exp(n,rlogis(n,...))

I'm not sure about the PERT beta -- a few seconds of web browsing
suggests that it's a shifted, scaled version of the beta.  Hence

a + (b-a)*rbeta(n,shape1,shape2)

  cheers
Ben Bolker

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


[R] Questions about sort data

2006-09-04 Thread Xiao Zhao
Dear R users,
I am doing my project which I want to plot a piecewise function, I knew
that I can use the command segments to plot. But the problem is I want
to use my real data which needs me to sort of my data by using the 'if
else'command, I use it 
If(t[i]36) lambda-0.5
Else lambda-0.2
The funny thing is when I look at my data set, it did not follow my
command to sort data, also the final numbers do not change either.
The other problem is I have a big data set, sometimes I want to know the
number of some data, such as how many numbers are over 60, how to do
this by using R?
Does anybody got any ideas about my two pros?

Thank you in advance
Best wishes
Nessie^_^

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


Re: [R] Questions about sort data

2006-09-04 Thread Gabor Grothendieck
Regarding your other questions
On 9/4/06, Xiao Zhao [EMAIL PROTECTED] wrote:
 Dear R users,
 I am doing my project which I want to plot a piecewise function, I knew
 that I can use the command segments to plot. But the problem is I want
 to use my real data which needs me to sort of my data by using the 'if
 else'command, I use it
 If(t[i]36) lambda-0.5
 Else lambda-0.2
 The funny thing is when I look at my data set, it did not follow my
 command to sort data, also the final numbers do not change either.
 The other problem is I have a big data set, sometimes I want to know the
 number of some data, such as how many numbers are over 60, how to do
 this by using R?
 Does anybody got any ideas about my two pros?


Using the built in data set rivers there are

  sum(rivers  1000)

rivers greater than 1000 miles long.

See the last two lines on every message to r-help.

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


[R] Problems with 2-D Kernel Density Estimation using MASS KernSmooth

2006-09-04 Thread Chanchal Kumar
Dear R-Users,
   
I am using the two dimensional Kernel density estimation function in
MASS package (specifically kde2d) and am having a recurrent problem.
The problem is that when my data vectors have less then 3 entries
then the functions gives me the density estimate. But when the vector
size increases beyond 3 then I get an error. I have pasted below my
steps and the error message. 

Pmin-2
Pmax-14
Mmin-100
Mmax-1000
N-200

dens-kde2d(data3[[3]],data3[[4]], h = c(width.SJ(data3[[3]],nb=100,
+method=dpi), width.SJ(data3[[4]], nb=100,method=dpi)), n=N,
+lims=c(Pmin,Pmax,log10(Mmin),log10(Mmax)))

When data3[[3]] and data3[[4]] are each less then 3 entries then the
function runs fine else it gives me following error message:

Error in SDh(cnt, (2.394/(n * TD))^(1/7), n, d) : NA/NaN/Inf in foreign
function call (arg 5)

I also tried using the bkde2D from KernSmooth package and am getting
the same error for longer vector sizes. 

I shall be thankful if you could suggest what is going wrong.
Alternatively I will be glad to know about other possible ways of fast
density estimation on bigger vectors. Thanks in advance!

Best Regards,
Chanchal 
===
Chanchal Kumar, Ph.D. Candidate
Dept. of Proteomics and Signal Transduction
Max Planck Institute of Biochemistry
Am Klopferspitz 18
82152 D-Martinsried (near Munich)
Germany
e-mail: [EMAIL PROTECTED]
Phone: (Office) +49 (0) 89 8578 2296
Fax:(Office) +49 (0) 89 8578 2219
http://www.biochem.mpg.de/mann/

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


Re: [R] Cross-correlation between two time series data

2006-09-04 Thread Spencer Graves
  'ccf' provides a plot with red, dashed lines indicating an 
approximate 95% threshold for the correlation. 

  Beyond that, with any particular model fit, you can get confidence 
intervals and anova tests for any particular parameter estimated. 

  If neither of these are adequate, I suppose one might be able to 
try Markov Chain Monte Carlo, but I've never used that, so I can't 
comment further on that. 

  If you would like more help from this listserve, please provide 
more detail of your application including commented, minimal, 
self-contained, reproducible code, explaining something you've tried and 
why it is not adequate (as suggested in the posting guide 
www.R-project.org/posting-guide.html). 

  Hope this helps. 
  Spencer Graves

Juni Joshi wrote:
Hi all,

I  have  two  time  series  data  (say  x  and  y). I am interested to
calculate the correlation between them and its confidence interval (or
to  test  no  correlation). Function cor.test(x,y) does the test of no
correlation. But this test probably is wrong because of autocorrelated
data.

ccf()  calculates the correlation between two series data. But it does
not  provide  the  confidence intervals of cross correlation. Is there
any  function  that  calculates the confidence interval of correlation
between  two  time  series data or performs the test of no correlation
between two time series data.

Thanks.

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


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


[R] opening files in directory

2006-09-04 Thread Ffenics
Hi there
I want to be able to take all the files in a given directory, read them in one 
at a time, calculate a distance matrix for them (the files are data matrices) 
and then print them out to separate files. This is the code I thought I would 
be able to use
(all files are in directory data_files)
for(i in 1:length(files))
+ {
+ x-read.table(data_files/files[[i]])
+ dist-dist(x, method=euclidean, diag=TRUE)
+ mat-as.matrix(dist)
+ write.table(mat, file=files[[i]])
+ }
But I get this error when I try to open the first file using read.table
Error in file(file, r) : unable to open connection
In addition: Warning message:
cannot open file 'data_files/files[[i]]'
if I try the read.table command without the quotation marks like so
x-read.table(data_matrix_files/files[[i]])
I get the error 
Error in read.table(data_matrix_files/files[[i]]) :
Object data_matrix_files not found
But if I go to the directory where the files are kept before starting up R, the 
read.table command without the quotation marks works.
I don't want to start up R in the same directory as the where the files I will 
be using reside though so how do I rectify this?
Any help much appreciated

[[alternative HTML version deleted]]

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


Re: [R] opening files in directory

2006-09-04 Thread Uwe Ligges
FAQ ...

Uwe Ligges



Ffenics wrote:
 Hi there
 I want to be able to take all the files in a given directory, read them in 
 one at a time, calculate a distance matrix for them (the files are data 
 matrices) and then print them out to separate files. This is the code I 
 thought I would be able to use
 (all files are in directory data_files)
 for(i in 1:length(files))
 + {
 + x-read.table(data_files/files[[i]])
 + dist-dist(x, method=euclidean, diag=TRUE)
 + mat-as.matrix(dist)
 + write.table(mat, file=files[[i]])
 + }
 But I get this error when I try to open the first file using read.table
 Error in file(file, r) : unable to open connection
 In addition: Warning message:
 cannot open file 'data_files/files[[i]]'
 if I try the read.table command without the quotation marks like so
 x-read.table(data_matrix_files/files[[i]])
 I get the error 
 Error in read.table(data_matrix_files/files[[i]]) :
 Object data_matrix_files not found
 But if I go to the directory where the files are kept before starting up R, 
 the read.table command without the quotation marks works.
 I don't want to start up R in the same directory as the where the files I 
 will be using reside though so how do I rectify this?
 Any help much appreciated
 
   [[alternative HTML version deleted]]
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

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


[R] read csv

2006-09-04 Thread Georg Otto

Hi,

I have a csv file where the number of filled columns varies in the
different rows:

Sun 5-Feb-06,15,,,01:30:00,0:06:00,
Mon 6-Feb-06,,
Tue 7-Feb-06,7,,,00:41:00,0:05:51,
Wed 8-Feb-06,,

I would like to use read.table (or whatever is appropriate) to read in
only those rows that have two or more columns filled. Any hint will be
appreciated.

Georg

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


Re: [R] opening files in directory

2006-09-04 Thread Ffenics
Yes, already looked at the FAQ. Thats how I have got where I am.

I didnt see a solution to this particular problem on there though

Thanks anyway.

Uwe Ligges [EMAIL PROTECTED] wrote: FAQ ...

Uwe Ligges



Ffenics wrote:
 Hi there
 I want to be able to take all the files in a given directory, read them in 
 one at a time, calculate a distance matrix for them (the files are data 
 matrices) and then print them out to separate files. This is the code I 
 thought I would be able to use
 (all files are in directory data_files)
 for(i in 1:length(files))
 + {
 + x-read.table(data_files/files[[i]])
 + dist-dist(x, method=euclidean, diag=TRUE)
 + mat-as.matrix(dist)
 + write.table(mat, file=files[[i]])
 + }
  But I get this error when I try to open the first file using read.table
 Error in file(file, r) : unable to open connection
 In addition: Warning message:
 cannot open file 'data_files/files[[i]]'
 if I try the read.table command without the quotation marks like so
 x-read.table(data_matrix_files/files[[i]])
 I get the error 
 Error in read.table(data_matrix_files/files[[i]]) :
 Object data_matrix_files not found
 But if I go to the directory where the files are kept before starting up R, 
 the read.table command without the quotation marks works.
 I don't want to start up R in the same directory as the where the files I 
 will be using reside though so how do I rectify this?
 Any help much appreciated
 
  [[alternative HTML version deleted]]
 
 __
 R-help@stat.math.ethz.ch mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.



[[alternative HTML version deleted]]

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


Re: [R] read csv

2006-09-04 Thread Prof Brian Ripley
On Mon, 4 Sep 2006, Georg Otto wrote:

 
 Hi,
 
 I have a csv file where the number of filled columns varies in the
 different rows:
 
 Sun 5-Feb-06,15,,,01:30:00,0:06:00,
 Mon 6-Feb-06,,
 Tue 7-Feb-06,7,,,00:41:00,0:05:51,
 Wed 8-Feb-06,,
 
 I would like to use read.table (or whatever is appropriate) to read in
 only those rows that have two or more columns filled. Any hint will be
 appreciated.

See ?read.table is the main hint.

Use read.table(fill=TRUE) and post-process, e.g. by

A - A[rowSums(!is.na(A))  2), ]

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

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


Re: [R] read csv

2006-09-04 Thread Thomas Hunger
 I have a csv file where the number of filled columns
 varies in the different rows:

I would hack it like this, but then I am totally new to R, 
which means you should not trust me.:

d - read.csv(testdata, header=F)

selection - apply(d, c(1), 
  function(x) {sum(!is.na(x)  x != )  2})

as.data.frame(t(as.data.frame(t(d))[selection]))


Tom

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


Re: [R] Fitting Pareto distribution to some data

2006-09-04 Thread William Asquith
Paul,
Package lmomco fits generalized pareto (three parameter) using method  
of L-moments.  I suspect that other packages that Brian identified  
use method of moments or other.

William

On Sep 3, 2006, at 10:55 AM, Paul Smith wrote:

 Dear All

 I am trying to fit Pareto distribution to some data. MASS package does
 not support Pareto distribution. Is there some alternative way?

 Thanks in advance,

 Paul

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

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


Re: [R] Subsetting vectors based on condition

2006-09-04 Thread J. Hosking
Mahesh Krishnan wrote:
 Hello,
 
 I have a question regarding subsetting of vectors. Here's an example of
 what I'm trying to do:
 
 vect.1 - c(76,195, 290, 380)
 
 vect.2 -  c(63,  95, 133, 170, 215, 253, 285, 299, 325, 375)
 
 I would like to subset vect.2 so that it has the same length as vect.1,
 and its numbers are the first corresponging higher value compared to
 vect.1.
 
 The output should be:
 
  final.output = (95, 215, 299, NA)
 
 What is the fastest/most eficient way to accompllish this in R?
 
 Thanks for the help in advance,
 
 Mahesh Krishnan
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 

  vect.1 - c(76,195, 290, 380)
  vect.2 -  c(63,  95, 133, 170, 215, 253, 285, 299, 325, 375)
  vect.2[ findInterval(vect.1,vect.2) + 1 ]
[1]  95 215 299  NA

J. R. M. Hosking

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


Re: [R] Fitting Pareto distribution to some data

2006-09-04 Thread Paul Smith
On 9/4/06, William Asquith [EMAIL PROTECTED] wrote:
 Package lmomco fits generalized pareto (three parameter) using method
 of L-moments.  I suspect that other packages that Brian identified
 use method of moments or other.

That is excellent to learn that, William. Thanks.

Paul


 On Sep 3, 2006, at 10:55 AM, Paul Smith wrote:

  Dear All
 
  I am trying to fit Pareto distribution to some data. MASS package does
  not support Pareto distribution. Is there some alternative way?
 
  Thanks in advance,
 
  Paul
 
  __
  R-help@stat.math.ethz.ch mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide http://www.R-project.org/posting-
  guide.html
  and provide commented, minimal, self-contained, reproducible code.



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


Re: [R] Optimization

2006-09-04 Thread Spencer Graves
  Have you considered talking logarithms of the expression you 
mentioned: 

  log(Yield) = a1*log(A)+b1*log(B)+c2*log(C)+...

where a1 = a/(a+b+...), etc.  This model has two constraints not present 
in ordinary least squares:  First, the intercept is assumed to be zero.  
Second, the coefficients in this log formulation must sum to 1.  If I 
were you, I might use something like lm to test them both. 

  To explain how, I'll modify the notation, replacing A by X1, B by 
X2, ..., up to Xkm1 (= X[k-1]) and Xk for k different environmental 
variables.  Then I might try something like the following: 

  fit0 - lm(log(Yield) ~ log(X1) + ... + log(Xk)-1 )
  fit1 - lm(log(Yield) ~ log(X1) + ... + log(Xk) )
  fit.1 - lm(log(Yield/Xk) ~ log(X1/Xk) + ... + log(Xkm1/Xk) )
  fit.0 - lm(log(Yield/Xk) ~ log(X1/Xk) + ... + log(Xkm1/Xk)-1 )

  anova(fit1, fit0) would test the no-constant model, and if I 
haven't made a mistake in this, anova(fit0, fit.0) and anova(fit1, 
fit.1) would test the constraint that all the coefficients should sum to 
1. 

  If you would like further help from this listserve, please provide 
commented, minimal, self-contained, reproducible code to help potential 
respondents understand your question and concerns (as suggested in the 
posting guide www.R-project.org/posting-guide.html). 

  Hope this helps. 
  Spencer Graves

Simone Vincenzi wrote:
 Dear R-list,
 I'm trying to estimate the relative importance of 6 environmental variables
 in determining clam yield. To estimate clam yield a previous work used the
 function Yield = (A^a*B^b*C^c...)^1/(a+b+c+...) where A,B,C... are the
 values of the environmental variables and the weights a,b,c... have not been
 calibrated on data but taken from literature. Now I'd like to estimate the
 weights a,b,c... by using a dataset with 110 observations of yield and
 values of the environmental variables. I'm wondering if it is feasible or if
 the number of observation is too low, if some data transformation is needed
 and which R function is the most appropriate to try to estimate the weights.
 Any help would be greatly appreciated.

 Simone Vincenzi 

 _
 Simone Vincenzi, PhD Student 
 Department of Environmental Sciences
 University of Parma
 Parco Area delle Scienze, 33/A, 43100 Parma, Italy
 Phone: +39 0521 905696
 Fax: +39 0521 906611
 e.mail: [EMAIL PROTECTED] 



 --

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


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


Re: [R] opening files in directory

2006-09-04 Thread Mike Nielsen
R won't do variable interpolation inside quotation marks as perl does.

You could try amending your code with, for e.g.

file.name-paste(sep=/,data_files,files[[i]])
x-read.table(file.name)

Regards,

Mike

On 9/4/06, Ffenics [EMAIL PROTECTED] wrote:
 Hi there
 I want to be able to take all the files in a given directory, read them in 
 one at a time, calculate a distance matrix for them (the files are data 
 matrices) and then print them out to separate files. This is the code I 
 thought I would be able to use
 (all files are in directory data_files)
 for(i in 1:length(files))
 + {
 + x-read.table(data_files/files[[i]])
 + dist-dist(x, method=euclidean, diag=TRUE)
 + mat-as.matrix(dist)
 + write.table(mat, file=files[[i]])
 + }
 But I get this error when I try to open the first file using read.table
 Error in file(file, r) : unable to open connection
 In addition: Warning message:
 cannot open file 'data_files/files[[i]]'
 if I try the read.table command without the quotation marks like so
 x-read.table(data_matrix_files/files[[i]])
 I get the error
 Error in read.table(data_matrix_files/files[[i]]) :
 Object data_matrix_files not found
 But if I go to the directory where the files are kept before starting up R, 
 the read.table command without the quotation marks works.
 I don't want to start up R in the same directory as the where the files I 
 will be using reside though so how do I rectify this?
 Any help much appreciated

 [[alternative HTML version deleted]]

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




-- 
Regards,

Mike Nielsen

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


Re: [R] opening files in directory

2006-09-04 Thread Gabor Grothendieck
On 9/4/06, Mike Nielsen [EMAIL PROTECTED] wrote:
 R won't do variable interpolation inside quotation marks as perl does.

Just as an aside, gsubfn in package gsubfn will do perl-style
(well, sort of) string interpolation:

 library(gsubfn)
 i - 1
 gsubfn(x = data_files/file$i)
[1] data_files/file1

cati and cati0 in the same package provide for such interpolation
within cat.

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


Re: [R] abline and plot(augPred) help

2006-09-04 Thread Paul Murrell
Hi


Petr Pikal wrote:
 Dear all
 
 as I did not get any response on my post about abline and 
 plot(augPred)) I try again. I hope I do not break some posting guide 
 rules. I would try to contact package maintainer directly but there 
 is stated to be R-core people, so I feel R-help list shall be OK.
 
 I need to draw straight lines through augPred plotted panels 
 (vertical or horizontal) at specified point. I know I shall probably 
 use panel.abline but I am missing correct syntax. Below you can see 
 my attempts together with results. I hope somebody can point me to 
 right direction.
 
 I am probably somewhere close but I have no clue, which parameter I
 shall modify to get measured points, fitted lines and vertical lines
 in panels together.


The problem is that you do not know about the default panel function
that nlme:::plot.augPred() uses, so your panel functions are not
replicating all of the default behaviour as well as adding your vertical
lines.  Some possible solutons suggested below ...


 fm1 - lme(Orthodont)
 
 # standard plot
 plot(augPred(fm1, level = 0:1, length.out = 2))
 
 #plot with vertical but without points and fitted lines
 plot(augPred(fm1, level = 0:1, length.out = 2),
 panel=function(v,...) {
 panel.abline(v=10)}
 )
 
 # plot with vertical but without fitted lines
 plot(augPred(fm1, level = 0:1, length.out=2),
 panel=function(x,y,...) {
 panel.xyplot(x,y,...)
 panel.abline(v=10)}
 )
 
 # plot with vertical and with all points (fitted lines are drawn as 
 points)
 plot(augPred(fm1, level = 0:1),
 panel=function(x,y,...) {
 panel.xyplot(x,y,...)
 panel.abline(v=10)}
 )

One option is to take a sneak a peek at nlme:::plot.augPred() to see
what the default panel function is doing.  Here I have replicated the
default panel function and added a call to panel.abline().

plot(augPred(fm1, level = 0:1, length.out = 2),
  panel=function(x, y, subscripts, groups, ...) {
orig - groups[subscripts] == original
panel.xyplot(x[orig], y[orig], ...)
panel.superpose(x[!orig], y[!orig], subscripts[!orig],
groups, ..., type = l)
panel.abline(v=10)
  })

The problem with this approach is that you need to crawl around in the
code of nlme:::plot.augPred().  An alternative approach is to annotate
the plot after-the-fact.  This is shown below.

plot(augPred(fm1, level = 0:1, length.out = 2))
for (i in 1:5) {
  for (j in 1:6) {
if (i  5 || j  4) {
  trellis.focus(panel, j, i, highlight=FALSE)
  panel.abline(v=10)
}
  }
}

This avoids crawling around in code, but the problem with this is
knowing how many rows and columns of panels there are.  If you
explicitly controlled the 'layout' of the original plot, you could
guarantee that your annotation works properly.

Hope that helps.

Paul
-- 
Dr Paul Murrell
Department of Statistics
The University of Auckland
Private Bag 92019
Auckland
New Zealand
64 9 3737599 x85392
[EMAIL PROTECTED]
http://www.stat.auckland.ac.nz/~paul/

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


Re: [R] abline and plot(augPred) help

2006-09-04 Thread Gabor Grothendieck
On 9/4/06, Paul Murrell [EMAIL PROTECTED] wrote:
 Hi


 Petr Pikal wrote:
  Dear all
 
  as I did not get any response on my post about abline and
  plot(augPred)) I try again. I hope I do not break some posting guide
  rules. I would try to contact package maintainer directly but there
  is stated to be R-core people, so I feel R-help list shall be OK.
 
  I need to draw straight lines through augPred plotted panels
  (vertical or horizontal) at specified point. I know I shall probably
  use panel.abline but I am missing correct syntax. Below you can see
  my attempts together with results. I hope somebody can point me to
  right direction.
 
  I am probably somewhere close but I have no clue, which parameter I
  shall modify to get measured points, fitted lines and vertical lines
  in panels together.


 The problem is that you do not know about the default panel function
 that nlme:::plot.augPred() uses, so your panel functions are not
 replicating all of the default behaviour as well as adding your vertical
 lines.  Some possible solutons suggested below ...


  fm1 - lme(Orthodont)
 
  # standard plot
  plot(augPred(fm1, level = 0:1, length.out = 2))
 
  #plot with vertical but without points and fitted lines
  plot(augPred(fm1, level = 0:1, length.out = 2),
  panel=function(v,...) {
  panel.abline(v=10)}
  )
 
  # plot with vertical but without fitted lines
  plot(augPred(fm1, level = 0:1, length.out=2),
  panel=function(x,y,...) {
  panel.xyplot(x,y,...)
  panel.abline(v=10)}
  )
 
  # plot with vertical and with all points (fitted lines are drawn as
  points)
  plot(augPred(fm1, level = 0:1),
  panel=function(x,y,...) {
  panel.xyplot(x,y,...)
  panel.abline(v=10)}
  )

 One option is to take a sneak a peek at nlme:::plot.augPred() to see
 what the default panel function is doing.  Here I have replicated the
 default panel function and added a call to panel.abline().

 plot(augPred(fm1, level = 0:1, length.out = 2),
  panel=function(x, y, subscripts, groups, ...) {
orig - groups[subscripts] == original
panel.xyplot(x[orig], y[orig], ...)
panel.superpose(x[!orig], y[!orig], subscripts[!orig],
groups, ..., type = l)
panel.abline(v=10)
  })

 The problem with this approach is that you need to crawl around in the
 code of nlme:::plot.augPred().  An alternative approach is to annotate
 the plot after-the-fact.  This is shown below.

 plot(augPred(fm1, level = 0:1, length.out = 2))
 for (i in 1:5) {
  for (j in 1:6) {
if (i  5 || j  4) {
  trellis.focus(panel, j, i, highlight=FALSE)
  panel.abline(v=10)
}
  }
 }

 This avoids crawling around in code, but the problem with this is
 knowing how many rows and columns of panels there are.  If you
 explicitly controlled the 'layout' of the original plot, you could
 guarantee that your annotation works properly.


You can find that out with trellis.currentLayout:

tcL - trellis.currentLayout()
for(i in 1:nrow(tcL))
  for(j in 1:ncol(tcL))
if (tcL[i,j]  0) {
trellis.focus(panel, j, i, highlight = FALSE)
panel.abline(v = 10)
trellis.unfocus()
}

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


[R] Coercing elements of a matrix from integer to double

2006-09-04 Thread Thaden, John J
Ive been converting elements of matrices and arrays, e.g.,  from
Integers to double-precision, by vectorizing the matrix and then
remaking it. Alternatively, I can redefine one element as double 
which then redefines them all.  Both methods are quick, so I guess
I shouldn't complain, but I would have thought there'd be something
more obvious. Have I missed it?

Here's my redimensioning example:  

## Matrix M...
M - 1:2e6 ; dim(Mi) - c(1e3,2e3)
dim(M)
class(M)
## ...has integer elements, e.g.,
class(M[1,1])

## The as.double() command changes
## these to double-precision, but it
## also strips away dimensions...
Md - as.double(Mi)
dim(Md)
class(Md)
## ...so I have to put them back.
dim(Md) - dim(Mi)
dim(Md) 
class(Md)
class(Md[1,1])

Here's my tail wagging the dog example:

M[1,1] - as.double(M[1,1])
class(M[2,2])

Thanks,
-John Thaden

Confidentiality Notice: This e-mail message, including any a...{{dropped}}

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


Re: [R] Coercing elements of a matrix from integer to double

2006-09-04 Thread Peter Dalgaard
Thaden, John J [EMAIL PROTECTED] writes:

 Ive been converting elements of matrices and arrays, e.g.,  from
 Integers to double-precision, by vectorizing the matrix and then
 remaking it. Alternatively, I can redefine one element as double 
 which then redefines them all.  Both methods are quick, so I guess
 I shouldn't complain, but I would have thought there'd be something
 more obvious. Have I missed it?

storage.mode(M) - double


-- 
   O__   Peter Dalgaard Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark  Ph:  (+45) 35327918
~~ - ([EMAIL PROTECTED])  FAX: (+45) 35327907

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


Re: [R] how to fit gauss beam?

2006-09-04 Thread Thomas Hunger
Reply to self:
 I am having a hard time fitting a gauss beam using R. In
 gnutplot I did something like

 $ w(z) = w0 * sqrt(1+(z/z0)**2)
 $ fit w(z) 'before_eom.txt' using 1:2 via w0, z0

 to obtain w0 and z0. Now I want to do the same in R. I
 tried a linear model like this (r = radius, z =

This works fine:

data - read.table (nach_eom.tab, header=T)

M - 1.1
nls (major ~ M*w0*sqrt(1+(d/z0)^2), 
 data = data, 
 start = list(w0 = 460, z0=0.02))

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


Re: [R] Cross-correlation between two time series data

2006-09-04 Thread Andrew Robinson
Jun,

If your interest is to estimate the correlation and either a
confidence interval or a test for no correlation, then you might try
to proceed as follows.  This is a Monte-Carlo significance test, and a
useful strategy.

1) use ccf() to compute the cross-correlation between x and y.

2) repeat the following steps, say, 1000 times.

2a) randomly reorder the values of one of the time series, say x.
Call the randomly reordered series x'. 

2b) use ccf() to compute the cross-correlation between x' and y.
Store that cross-correlation.

3) the 1000 cross-correlation estimates computed in step 2 are all
   estimating cross-correlation 0, conditional on the data.  A
   two-tailed test then is: if the cross-correlation computed in step
   1 is outside the (0.025, 0.975) quantiles of the empirical
   distribution of the cross-correlations computed in step 2, then,
   reject the null hypothesis that x and y are uncorrelated, with size
   0.05.

I hope that this helps.

Andrew


Juni Joshi wrote:
Hi all,

I  have  two  time  series  data  (say  x  and  y). I am interested to
calculate the correlation between them and its confidence interval (or
to  test  no  correlation). Function cor.test(x,y) does the test of no
correlation. But this test probably is wrong because of autocorrelated
data.

ccf()  calculates the correlation between two series data. But it does
not  provide  the  confidence intervals of cross correlation. Is there
any  function  that  calculates the confidence interval of correlation
between  two  time  series data or performs the test of no correlation
between two time series data.

Thanks.

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


-- 
Andrew Robinson  
Department of Mathematics and StatisticsTel: +61-3-8344-9763
University of Melbourne, VIC 3010 Australia Fax: +61-3-8344-4599
Email: [EMAIL PROTECTED] http://www.ms.unimelb.edu.au

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


Re: [R] Subsetting vectors based on condition

2006-09-04 Thread Mahesh Krishnan


Sincere thanks to Jim Holtman and J. Hosking for their suggestions.
both their solutions work perfectly,
in particular the findInterval function is  what I was looking for.

Cheers.

Mahesh Krishnan

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


[R] Fitting generalized additive models with constraints?

2006-09-04 Thread David Reiss
Hello,
I am trying to fit a GAM for a simple model, a simple model, y ~ s(x0) +
s(x1) ; with a constraint that the fitted smooth functions s(x0) and s(x1)
have to each always be 0.

From the library documentation and a search of the R-site and R-help
archives I have not been able to decipher whether the following is possible
using this, or other GAM libraries, or whether I will have to try to roll
my own. I see from the mgcv docs that GAMs need to be constrained such that
the smooth functions have zero mean. Is there a way around this?

Is such a constraint possible?
thanks very much for any advice or pointers.
-David

[[alternative HTML version deleted]]

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


[R] Problem with Variance Components (and general glmm confusion)

2006-09-04 Thread Toby Gardner
Dear list, 

I am having some problems with extracting Variance Components from a 
random-effects model: 

I am running a simple random-effects model using lme: 

model-lme(y~1,random=~1|groupA/groupB)

which returns the output for the StdDev of the Random effects, and model AIC 
etc as expected. 

Until yesterday I was using R v. 2.0, and had no problem in calling the 
variance components of the above model using VarCorr(model), together with 
their 95% confidence intervals using intervals() - although for some response 
variables a call to intervals() returns the error: Cannot get confidence 
intervals on var-cov components: Non-positive definite approximate 
variance-covariance.  

I have now installed R v. 2.3.1 and am now experiencing odd behaviour with 
VarCorr(lme.object), with an error message typically being returned: 

Error in VarCorr(model) : no direct or inherited method for function 'VarCorr' 
for this call

Is this known to happen? For instance could it be due to the subsequent loading 
of new packages? (lme4 for instance?).  

To get around this problem I have tried running the same model using lmer: 

model2-lmer(y~1 + (1|groupA) + (1|groupB))

Should this not produce the same model? The variance components are very 
similar but not identical, making me think that I am doing something wrong. I 
am also correct in thinking that intervals() does not work with lmer? I get: 
Error in intervals(model2) : no applicable method for intervals 

GLMM

I have a general application question - please excuse my ignorance, I am 
relatively new to this and trying to find a way through the maze.  In short I 
need to compile generalized linear mixed models both for (a) Poisson data and 
(b) binonial data incorporating a two nested random factors, and I need to be 
able to extract AIC values as I am taking an information-theoretic approach to 
model selection.  Prior to sending an email to the list I have spent quite a 
few days reading the background on a number of functions, all of which offer 
potential for this; glmmML, glmmPQL, lmer, and glmmADMB.  I can understand that 
glmmPQL is unsuitable because there is no way of knowing the maximised 
likelihood, but is there much difference between the remaining three options? I 
have seen simulation comparisons published on this list between glmmADMB and 
glmmPQL and lmer, but it seems these are before the latest release of lmer, and 
also they do not evaluate glmmML.  To a newcomer this myriad !
 of options is bewildering, can anyone offer advice as to the most robust 
approach? 

Many thanks for your time and patience, 

Toby Gardner 

School of Environmental Sciences
University of East Anglia
Norwich, NR4 7TJ
United Kingdom
Email: [EMAIL PROTECTED]
Website: www.uea.ac.uk/~e387495

[[alternative HTML version deleted]]

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


Re: [R] Cross-correlation between two time series data

2006-09-04 Thread Spencer Graves
Hi, Andrew: 

  This will produce a permutation distribution for the correlation 
under the null hypothesis of zero correlation between the variables.  
This is a reasonable thing to do, and would probably produce limits more 
accurate than the dashed red lines on the 'ccf' plot.  However, they 
would NOT be confidence interval(s). 

  For a confidence interval on cross correlation, you'd have to 
hypothesize some cross correlation pattern between x and y, preferably 
parameterized parsimoniously, then somehow determine an appropriate 
range of values consistent with the data.  By the time you've done all 
that, you've effectively fit some model and constructed confidence 
intervals on the parameter(s). 

  Best Wishes,
  Spencer

Andrew Robinson wrote:
 Jun,

 If your interest is to estimate the correlation and either a
 confidence interval or a test for no correlation, then you might try
 to proceed as follows.  This is a Monte-Carlo significance test, and a
 useful strategy.

 1) use ccf() to compute the cross-correlation between x and y.

 2) repeat the following steps, say, 1000 times.

 2a) randomly reorder the values of one of the time series, say x.
 Call the randomly reordered series x'. 

 2b) use ccf() to compute the cross-correlation between x' and y.
 Store that cross-correlation.

 3) the 1000 cross-correlation estimates computed in step 2 are all
estimating cross-correlation 0, conditional on the data.  A
two-tailed test then is: if the cross-correlation computed in step
1 is outside the (0.025, 0.975) quantiles of the empirical
distribution of the cross-correlations computed in step 2, then,
reject the null hypothesis that x and y are uncorrelated, with size
0.05.

 I hope that this helps.

 Andrew


 Juni Joshi wrote:
   
Hi all,

I  have  two  time  series  data  (say  x  and  y). I am interested to
calculate the correlation between them and its confidence interval (or
to  test  no  correlation). Function cor.test(x,y) does the test of no
correlation. But this test probably is wrong because of autocorrelated
data.

ccf()  calculates the correlation between two series data. But it does
not  provide  the  confidence intervals of cross correlation. Is there
any  function  that  calculates the confidence interval of correlation
between  two  time  series data or performs the test of no correlation
between two time series data.

Thanks.

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

 



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


[R] plot a new picture against an old one to see the difference between them

2006-09-04 Thread Am Stat
Hello, useR:,

Suppose I have two plots made by using contour() function, say Cont1 and Cont2 
respectively. 

They have slightly difference because of the two slightly different data  I 
used.

I want to see the difference between them so I want to plot Cont2 on Cont1, are 
there any methods to plot it without filling the frame of Cont1 totally of 
Cont2.
I mean, how I can integreate the two plots together that they kind of have 
weighted colors?

Thanks very much in Advance!

Leon
[[alternative HTML version deleted]]

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


Re: [R] plot a new picture against an old one to see the difference between them

2006-09-04 Thread Am Stat
Dear Roger,

Thanks, that's really helpful,  do you know how to deal with it if the two 
plots are  generated by plot(), not by contour().

Best,

Leon


- Original Message - 
From: roger koenker [EMAIL PROTECTED]
To: Am Stat [EMAIL PROTECTED]
Sent: Monday, September 04, 2006 8:06 PM
Subject: Re: [R] plot a new picture against an old one to see the difference 
between them


 for the second call to contour  use the argument add=TRUE.

 On Sep 4, 2006, at 6:42 PM, Am Stat wrote:

 Hello, useR:,

 Suppose I have two plots made by using contour() function, say  Cont1 and 
 Cont2 respectively.

 They have slightly difference because of the two slightly different  data 
 I used.

 I want to see the difference between them so I want to plot Cont2  on 
 Cont1, are there any methods to plot it without filling the  frame of 
 Cont1 totally of Cont2.
 I mean, how I can integreate the two plots together that they kind  of 
 have weighted colors?

 Thanks very much in Advance!

 Leon
 [[alternative HTML version deleted]]

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


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


[R] Help with plotmath

2006-09-04 Thread Kenneth Cabrera
Hi R users:

How can I have several subscript number with a comma in a plot.

I would like to have the LaTeX equivalent of

x_{i,j}.

I try:

plot(1:10,1:10,type=n)
text(5,5,expression(x[i,j]))

but it doesn´t work.

Thank you for your help.

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


[R] Sweave and the [ function

2006-09-04 Thread Ross Darnell
I am wanting to use the [ operator in an S-chunk, e.g.


=
str(women)
women$height
women[,1]
[(women,1)
@

to show the equivalence of  three methods of extracting an element from 
a data.frame.

However Sweave returns the last of these as

women[1]

in the S input chunk

How can I force it not to do this and return [(women,1)

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


Re: [R] Question on Chi-square of null model in sem package

2006-09-04 Thread Guo Wei-Wei
Dear Pref. Fox

Sorry, I didn't receive your reply. I try the new sem package. It's
great. The following is the results that I got. The fit indices are
fine.

 Model Chisquare =  208   Df =  98 Pr(Chisq) = 6.6e-10
 Chisquare (null model) =  1741   Df =  120
 Goodness-of-fit index =  0.9
 Adjusted goodness-of-fit index =  0.87
 RMSEA index =  0.066   90 % CI: (0.054, 0.079)
 Bentler-Bonnett NFI =  0.88
 Tucker-Lewis NNFI =  0.92
 Bentler CFI =  0.93
 BIC =  -336

Thank you very much. You help me out so many problems.

Best wishes,
Wei-Wei


2006/9/4, John Fox [EMAIL PROTECTED]:
 Dear Wei-Wei,

 As I explained to you in private email yesterday (perhaps you didn't receive
 my reply?), the problem that you point out is due to a bug in the sem
 function that I fixed some time ago and then inadvertently reintroduced.
 Yesterday, I sent a corrected version of the sem package (0.9-5) to CRAN;
 the source package is there now and I'm sure that the compiled Windows
 package will appear in due course.

 Thank you once more for bringing the problem to my attention.

 John


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


Re: [R] Cross-correlation between two time series data

2006-09-04 Thread Andrew Robinson
Hi Spencer,

you are quite right.  I should have been careful to emphasize that the
strategy I suggested was intended only to produce the test for no
correlation clause of the either a confidence interval or a test for
no correlation sentence.

Cheers

Andrew

On Mon, Sep 04, 2006 at 04:00:44PM -0700, Spencer Graves wrote:
 Hi, Andrew: 
 
  This will produce a permutation distribution for the correlation 
 under the null hypothesis of zero correlation between the variables.  
 This is a reasonable thing to do, and would probably produce limits more 
 accurate than the dashed red lines on the 'ccf' plot.  However, they 
 would NOT be confidence interval(s). 
 
  For a confidence interval on cross correlation, you'd have to 
 hypothesize some cross correlation pattern between x and y, preferably 
 parameterized parsimoniously, then somehow determine an appropriate 
 range of values consistent with the data.  By the time you've done all 
 that, you've effectively fit some model and constructed confidence 
 intervals on the parameter(s). 
 
  Best Wishes,
  Spencer
 
 Andrew Robinson wrote:
 Jun,
 
 If your interest is to estimate the correlation and either a
 confidence interval or a test for no correlation, then you might try
 to proceed as follows.  This is a Monte-Carlo significance test, and a
 useful strategy.
 
 1) use ccf() to compute the cross-correlation between x and y.
 
 2) repeat the following steps, say, 1000 times.
 
 2a) randomly reorder the values of one of the time series, say x.
 Call the randomly reordered series x'. 
 
 2b) use ccf() to compute the cross-correlation between x' and y.
 Store that cross-correlation.
 
 3) the 1000 cross-correlation estimates computed in step 2 are all
estimating cross-correlation 0, conditional on the data.  A
two-tailed test then is: if the cross-correlation computed in step
1 is outside the (0.025, 0.975) quantiles of the empirical
distribution of the cross-correlations computed in step 2, then,
reject the null hypothesis that x and y are uncorrelated, with size
0.05.
 
 I hope that this helps.
 
 Andrew
 
 
 Juni Joshi wrote:
   
Hi all,
 
I  have  two  time  series  data  (say  x  and  y). I am interested to
calculate the correlation between them and its confidence interval (or
to  test  no  correlation). Function cor.test(x,y) does the test of no
correlation. But this test probably is wrong because of autocorrelated
data.
 
ccf()  calculates the correlation between two series data. But it does
not  provide  the  confidence intervals of cross correlation. Is there
any  function  that  calculates the confidence interval of correlation
between  two  time  series data or performs the test of no correlation
between two time series data.
 
Thanks.
 
Jun
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide 
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 
 
   

-- 
Andrew Robinson  
Department of Mathematics and StatisticsTel: +61-3-8344-9763
University of Melbourne, VIC 3010 Australia Fax: +61-3-8344-4599
Email: [EMAIL PROTECTED] http://www.ms.unimelb.edu.au

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


Re: [R] Help with plotmath

2006-09-04 Thread Gabor Grothendieck
Try:

text(5,5,expression(x[i * , * j]))

On 9/4/06, Kenneth Cabrera [EMAIL PROTECTED] wrote:
 Hi R users:

 How can I have several subscript number with a comma in a plot.

 I would like to have the LaTeX equivalent of

 x_{i,j}.

 I try:

 plot(1:10,1:10,type=n)
 text(5,5,expression(x[i,j]))

 but it doesn´t work.

 Thank you for your help.

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


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


Re: [R] Sweave and the [ function

2006-09-04 Thread hadley wickham
 =
 str(women)
 women$height
 women[,1]
 [(women,1)
 @

 to show the equivalence of  three methods of extracting an element from
 a data.frame.

 However Sweave returns the last of these as

 women[1]

 in the S input chunk

 How can I force it not to do this and return [(women,1)

I don't think you can.  Sweave parses your R code and from then on
uses the internal R representation.  R normalises the parse tree in
certain ways (eg. strips comments, formats source code, and clearly
normalises some function calls).  Since sweave uses this, and not the
original text, I don't think there is anyway to get around this,
unless there is some trick during parsing.

(And don't forget women[[1]])

Hadley

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


Re: [R] Sweave and the [ function

2006-09-04 Thread Vincent Goulet
Le Mardi 5 Septembre 2006 0:03, hadley wickham a écrit :
  =
  str(women)
  women$height
  women[,1]
  [(women,1)
  @
 
  to show the equivalence of  three methods of extracting an element from
  a data.frame.
 
  However Sweave returns the last of these as
 
  women[1]
 
  in the S input chunk
 
  How can I force it not to do this and return [(women,1)

 I don't think you can.  Sweave parses your R code and from then on
 uses the internal R representation.  R normalises the parse tree in
 certain ways (eg. strips comments, formats source code, and clearly
 normalises some function calls).  Since sweave uses this, and not the
 original text, I don't think there is anyway to get around this,
 unless there is some trick during parsing.

 (And don't forget women[[1]])

 Hadley

So here's a workaround (untested):

echo=TRUE, eval=TRUE=
str(women)
women$height
women[,1]
@
echo=TRUE, eval=FALSE=
[(women,1)
@
echo=FALSE, eval=TRUE=
[(women,1)
@

I often end up doing similar things.

HTHVincent

-- 
  Vincent Goulet, Professeur agrégé
  École d'actuariat
  Université Laval, Québec 
  [EMAIL PROTECTED]   http://vgoulet.act.ulaval.ca

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


[R] colorRamp

2006-09-04 Thread Rick Bischoff
Hi,

I am using colorRamp in the following way.  I am *sure* there is a  
better way to do this, so if you'd be so kind to show me the true R way:

Step 0: Create a new variable, say, x, that maps some other  
continuous variable I have onto the [0,1] line.
Step 1:  Store the result from colorRamp (a function), into, say, test

  test - colorRamp(mypalette)

Step 2:  In my data frame, data

  data$colorTemp - test(data$x)

Step 3: Write a new function

bob - function(temp) { rgb(temp[1],temp[2],temp[3],maxColorValue=255) }

Step 4:

  for (i in 1:dim(data)[1]) data[i,color] - bob(data[i,  
colorTemp])

Step 5:

map(states, region=data$region, fill=T, col=data$color)

Thanks in advance!
Rick

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


[R] Quick question about lm()

2006-09-04 Thread Tong Wang
Hi, 
 Feel awkward to ask , but really couldn't find a answer anywhere,   How 
could I extract the R^2 and t-stat. from the 
result of lm()?
 Thanks a lot. 

best

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


Re: [R] Quick question about lm()

2006-09-04 Thread Christos Hatzis
Say,

my.lm - lm(y ~ x, data=my.data)

Then if you try:

names(summary(my.lm)) 

you will see the components of the summary.lm object.  The coefficients and
t-statistics can be extracted by

summary(my.lm)$coefficients

and similarly for the r-squared and other statistics provided in the summary
report.

-Christos
 

-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Tong Wang
Sent: Tuesday, September 05, 2006 1:35 AM
To: r-help@stat.math.ethz.ch
Subject: [R] Quick question about lm()

Hi, 
 Feel awkward to ask , but really couldn't find a answer anywhere,   How
could I extract the R^2 and t-stat. from the 
result of lm()?
 Thanks a lot. 

best

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

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