[R] triangular matrices input/output

2012-05-16 Thread casperyc
Hi, Is there any package that deals with triangular matrices? Say ways of inputting an upper (lower) triangular matrix? Or convert a vector of length 6 to an upper (lower) triangular matrix (by row/column)? Thanks! - ## PhD candidate in Statistics Big R Fan Big LEGO

Re: [R] problem with Gauss Hermite ( x and w )

2012-05-10 Thread casperyc
Hi, I know what complex number are, but I am not sure what you meant by that? ##CODES### 2.5^(-2.4) [1] 0.1109032 -2.5^(-2.4) [1] -0.1109032 ##CODES### works fine. Negative powers mean they take the reciprocal and as far as I am concerned, real^real

Re: [R] problem with Gauss Hermite ( x and w )

2012-05-10 Thread casperyc
Rui Barradas wrote real^real is not necessarily real. The most well known example is (-1)^0.5 = imaginary unit. Damn, can't believe it! It's a silly mistake! Now that something wonders me is that when applying the Gaussian Hermit, sum w f(x_i) What happens when f(x_i) does not

[R] problem with Gauss Hermite ( x and w )

2012-05-09 Thread casperyc
Hi all, I am using the 'gaussHermite' function from the 'pracma' library CODES ### library(pracma) cc=gaussHermite(10) cc$x^2 cc$x^5 cc$x^4 CODES ### as far so good. However, it does NOT work for any NON integer values, say CODES

[R] something weird in integration (pracma library)

2012-05-05 Thread casperyc
Hi, library(pracma) k=20 mu=4.5 casigma=17000 myint=function(j) { quadinf(function(x) (1/(1+exp(-x)))^j*(1-1/(1+exp(-x)))^(k-j)*dnorm(x,mu,casigma),-Inf,Inf) } sapply(0:k,myint)

[R] Numerical integration again

2012-04-18 Thread casperyc
Hi all, Here is an integration function require(pracma) # for 'quadinf' myint=function(j) { quadinf(function(x) (1/(1+exp(-x)))^j*(1-1/(1+exp(-x)))^(k-j)*dnorm(x,mu,casigma),-Inf,Inf) } in any optimization routine. It works fine most of the time but failed with some particular sets of

[R] how to avoid this : underflow occurred in 'lgammacor'

2012-04-04 Thread casperyc
Hi all, I am constructing a likelihood involving the following lbeta(j + a, k - j + b) where j,k are constants and a and b are parameters (0). While doing the optimization, the error sometimes occurs, In lbeta(j + a, k - j + b) : underflow occurred in 'lgammacor' Is there a way to avoid it?

[R] Create Model Object (setClass?setMethod?)

2012-04-03 Thread casperyc
Hi all, I have a self written likelihood as a model and functions to optimize and get fitted values, confidence intervals ect. I wonder if there is a way to define a 'class', or a 'model' (or a certain object)? so that I can use 'summary' to produce a summary like it does for a lm object. Also,

Re: [R] R numerical integration

2012-03-25 Thread casperyc
The quadinf command in library pracma still fails when mu=-2.986731 with sigma=53415.18. While Maple gives me an estimate of 0.5001701024. Maple: (for those who are interested) myf:=(mu,sigma)-

[R] Vectorize (scalar) function

2012-03-23 Thread casperyc
Hi all, myint=function(mu,sigma){ integrate(function(x) dnorm(x,mu,sigma)/(1+exp(-x)),-Inf,Inf)$value } x=seq(0,50,length=3000) x=x[-1] plot(x,myint(4,x)) # not working yet I think I have to 'Vectorize' it somehow? What's a scalar function? and a primitive function? Thanks. casper

[R] show and produce PDF file with pdf() and dev.off( ) in function

2012-03-23 Thread casperyc
Hi all, I know how to use pdf() and dev.off() to produce and save a graph. However, when I put them in a function say myplot(x=1:20){ pdf(xplot.pdf) plot(x) dev.off() } the function work. But is there a way show the graph in R as well as saving it to the workspace? Thanks. casper

Re: [R] Quicker way to apply values to a function

2012-03-23 Thread casperyc
Hi Petr, Thanks for confirming that the integral is bounded. I was thinking about the same thing. However, this requires that 'sigma' is positive. The actual problem occurred in my optimization routine, where i have set the parameter sigma=exp(para), where para is the logit of a uniform

[R] R numerical integration

2012-03-23 Thread casperyc
Hi all, Is there any other packages to do numerical integration other than the default 'integrate'? Basically, I am integrating: integrate(function(x) dnorm(x,mu,sigma)/(1+exp(-x)),-Inf,Inf)$value The integration is ok provided sigma is 0. However, when mu=-1.645074 and sigma=17535.26 It

[R] Quicker way to apply values to a function

2012-03-22 Thread casperyc
Hi all, myint=function(mu,sigma){ integrate(function(x) dnorm(x,mu,sigma)/(1+exp(-x)),-Inf,Inf)$value } mymu=seq(-3,3,length(1000)) mysigma=seq(0,1,length(500))[-1] k=1 v=c() for (j in 1:length(mymu)) { for (i in 1:length(mysigma)) { v[k]=myint(mymu[j],mysigma[i])

[R] R CMD Rd2pdf [options] files

2012-03-14 Thread casperyc
Hi all, Sorry if I seem a bit pissed because I am! Basically, I have written some script and I wanted to make a documentation for it using .Rd format. I followed the example here in

Re: [R] function input as variable name (deparse/quote/paste) ??

2012-03-11 Thread casperyc
Thank you everyone for your reply. Like I said in my original post, this is just a demonstrative example of my 'big' self written script. My 'big' function take several inputs, of which the first 1 is the dataset and returns a LIST variable 'res-list()' to the workspace with many information.

[R] function input as variable name (deparse/quote/paste) ??

2012-03-10 Thread casperyc
Hi all Say I have a function: myname=function(dat,x=5,y=6){ res-x+y-dat } for various input such as myname(dat1) myname(dat2) myname(dat3) myname(dat4) myname(dat5) how should I modify the 'res' line, to have new informative variable name correspondingly, such as dat1.res dat2.res

Re: [R] function input as variable name (deparse/quote/paste) ??

2012-03-10 Thread casperyc
Sorry if I wasn't stating what I really wanted or it was a bit confusing. Basically, there are MANY datasets to run suing the same function I have written a function to analyze it and returns a LIST of useful out put in the variable 'res' (to the workspace). I also created another script run.r

[R] how to quote factors in a function?

2011-09-11 Thread casperyc
### code ### x=sample(LETTERS[1:26],100,T) prob=function(y){ count=sum(x==y) total=length(x) count/total } ### end ### How do I quote the letters A,B,C,D ect, in order to make the prob function return the weights of the desired Letter? Thanks! --

Re: [R] Binary response GLM Question

2011-06-05 Thread casperyc
Hi Josh, Thank you for your reply. The reason I thank Y (0 and 1) here as p is because I think each observation is just a bernulli trial, so in this case the binomial n=1. And yet R still fits it (with the logit link) . I know the expression for the logit link, so I assumed I can take y here as

[R] Binary response GLM Question

2011-06-04 Thread casperyc
Hi all, I have a problem with binary response data in GLM fitting. The problem is that the y take only 1 or 0, and if I use logit link, it is the log of the odds ratio, which is p/(1-p). In my situation, think y is p, so sometimes the odds is 0, sometimes it is 1/0, which is (should be) undefine?

[R] what controls the limits of x axis in hclust

2011-03-07 Thread casperyc
Hi all, I followed the example in http://www.statmethods.net/advstats/cluster.html and apply it to one of my own dataset, I got this tiny problem with the boreders, http://r.789695.n4.nabble.com/file/n3340238/temp.png The red rectangle are 'outside' the plot, so I want to know how do I

Re: [R] how to see what's wrong with a self written function?

2010-12-30 Thread casperyc
Hi Petr Savicky, It is very interesting and a good way to plot, so that helps me to visualized the method, also easier to see what's wrong. The readline is new to me, and I should learn how to use it in the future. Thanks! casper -- View this message in context:

[R] how to see what's wrong with a self written function?

2010-12-21 Thread casperyc
Hi all, I am writing a simple function to implement regularfalsi (secant) method. ### regulafalsi=function(f,x0,x1){ x=c() x[1]=x1 i=1 while ( f(x[i])!=0 ) { i=i+1 if (i==2) {

[R] use of 'apply' for 'hist'

2010-12-18 Thread casperyc
Hi all, ## dof=c(1,2,4,8,16,32) Q5=matrix(rt(100,dof),100,6,T,dimnames=list(NULL,dof)) par(mfrow=c(2,6)) apply(Q5,2,hist) myf=function(x){ qqnorm(x);qqline(x) } apply(Q5,2,myf) ## These looks ok. However, I would

[R] how to use paste in function for expression?

2010-12-18 Thread casperyc
Hi all, ## dof=c(1,2,4,8,16,32) Q5=matrix(rt(100,dof),100,6,T,dimnames=list(NULL,dof)) par(mfcol=c(2,6)) for (i in 1:6){ dof2=dof[i] hist(Q5[,i],main=paste(t[,dof2,],sep=)) qqnorm(Q5[,i]) qqline(Q5[,i]) }

Re: [R] use of 'apply' for 'hist'

2010-12-18 Thread casperyc
Thanks Dieter. Casper -- View this message in context: http://r.789695.n4.nabble.com/use-of-apply-for-hist-tp3093811p3093938.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list

Re: [R] use of 'apply' for 'hist'

2010-12-18 Thread casperyc
HI Sarah, I will just use a loop then. I think my colnames are fine. Thanks! casper -- View this message in context: http://r.789695.n4.nabble.com/use-of-apply-for-hist-tp3093811p3093937.html Sent from the R help mailing list archive at Nabble.com.

[R] how to use expression as function arguements?

2010-12-18 Thread casperyc
Hi all, # integ=function(f,n){ # computes vector y = f(x) x=runif(1) y=f hatI=mean(y) hatI } # example of use integ(x^2+x+100,100) # it returns an error says no obj 'x' how do I

Re: [R] How to print colorful R output??

2010-12-12 Thread casperyc
Hi All, My aim is actually not that complicated as you guys understand. What I want is this, when I print it by clicking File-- Print... It gaves me a black white output. But what I want is 'red', for all the codes i typed in, 'blue', for the R output, just like the console. Thanks! (I

[R] How to print colorful R output??

2010-12-10 Thread casperyc
Hi All, I wonder if there is a way to print the R output with COLOR? Not the color plots, but the outputs in the console. Thank. casper -- View this message in context: http://r.789695.n4.nabble.com/How-to-print-colorful-R-output-tp3082750p3082750.html Sent from the R help mailing list

[R] how to add these axis label?

2010-12-08 Thread casperyc
Hi All, How do I add these axis labels? ### p=seq(0,1,length.out=500) p=p[-c(1,length(p))] g1=log(p/(1-p)) g2=qnorm(p) g3=log(-log(1-p)) g4=-log(-log(p)) plot(p,g1, 'n',ylim=c(-5,5),las=1, bty='n', xaxt='n',yaxt='n',

Re: [R] how to add these

2010-12-08 Thread casperyc
Thomas Stewart wrote: It isn't clear to me what you want to do. Do you want the axes to show? Do you want labels for the lines? Do you want a legend? What is your desired output? -tgs On Wed, Dec 8, 2010 at 2:42 PM, casperyc caspe...@hotmail.co.uk wrote: Hi All, How do I

Re: [R] how to add these

2010-12-08 Thread casperyc
thanks! problem solved! casper -- View this message in context: http://r.789695.n4.nabble.com/how-to-add-these-tp3078908p3079340.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list

[R] how to divide each column in a matrix by its colSums?

2010-11-28 Thread casperyc
Hi, I have a matrix, say m=matrix(c( 983,679,134, 383,416,84, 2892,2625,570 ),nrow=3 ) i can find its row/col sum by rowSums(m) colSums(m) How do I divide each row/column by its rowSum/colSums and still return in the matrix form? (i.e. the new rowSums/colSums

Re: [R] how to divide each column in a matrix by its colSums?

2010-11-28 Thread casperyc
In that case, there are values 1, which is clearly not what I wanted. Thanks. I think I should use prop.table -- View this message in context: http://r.789695.n4.nabble.com/how-to-divide-each-column-in-a-matrix-by-its-colSums-tp3062739p3062752.html Sent from the R help mailing list archive

Re: [R] how to divide each column in a matrix by its colSums?

2010-11-28 Thread casperyc
I am using prop.table(m,1) and prop.table(m,2) my aim. which I think is the most 'easy' way. Thanks. Casper -- View this message in context: http://r.789695.n4.nabble.com/how-to-divide-each-column-in-a-matrix-by-its-colSums-tp3062739p3062797.html Sent from the R help mailing list archive

[R] how to add frequencies to barplot

2010-11-20 Thread casperyc
Hi, I have count data x2=rep(c(0:3),c(13,80,60,27)) x2 0 1 2 3 13 80 60 27 I want to graph to be ploted as barplot(table(x2),density=4) how do I add relative frequency to it, like in hist(x2,labels=T) above the 'bar's Thanks. casper -- View this message in context:

Re: [R] how exactly does 'identify' work?

2010-11-18 Thread casperyc
Hi, I think the problem is 1 - when a linear model is fitted, ploting the qqnorm( test.lm$ res ) we dont 'know' what values are actually being used on the y-axis; and how do we refer to the ‘Index’ on the x-axis?? therefore, i dont know how to refer to the x and y coordinates in the

Re: [R] how exactly does 'identify' work?

2010-11-18 Thread casperyc
yes, i tried to modify the 2L part in plot.lm ### if (show[2L]) { ylim - range(rs, na.rm = TRUE) ylim[2L] - ylim[2L] + diff(ylim) * 0.075 qq - qqnorm(rs, main = main, ylab = ylab23, ylim = ylim, ...) if (qqline)

Re: [R] how exactly does 'identify' work?

2010-11-18 Thread casperyc
Omg! yes, it is working now! tmp - qqnorm( resid(test.lm) ) What a simple nice trick!!! Actually, i wasnt looking for the 'i'th label, I was looking for the 'row.names' as label, like I stated in the 1st post. identify(tmp, , row.names(test) ) is the label i have been trying to get.

[R] plot linear model problem

2010-11-16 Thread casperyc
Hi all, Say I fit a linear model, and saved it as 'test.lm' Then if I use plot(test.lm) it gives me 4 graphs How do I ask for a 'subset' of it?? say just want the 1st graph, the residual vs fitted values, or the 1,3,4th graph? I think I can use plot(test.lm[c(1,3,4)]) before, but now,

Re: [R] plot linear model problem

2010-11-16 Thread casperyc
Thank you both. casper -- View this message in context: http://r.789695.n4.nabble.com/plot-linear-model-problem-tp3045763p3045932.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list

[R] how exactly does 'identify' work?

2010-11-16 Thread casperyc
Hi all, # test=data.frame(x=1:26,y=-23.5+0.45*(1:26)+rnorm(26)) rownames(test)=LETTERS[1:26] attach(test) #test test.lm=lm(y~x) plot(test.lm,2) identify(test.lm$res,,row.names(test)) # not working plot(x,y) identify(x,y,row.names(test)) # works fine

[R] what's wrong with this 'length' in function?

2010-11-11 Thread casperyc
Hi all, I am having a trouble with this function I wrote ### p26=function(x,alpha){ # dummy variable j=1 ci=matrix(ncol=2,nrow=3) while (j4){ if (j==2) {x=x+c(-1,1)*0.5}

[R] par mfrow in function problem

2010-11-10 Thread casperyc
Hi all, I defined the following # myhist=function(x){ hist(x,xlab=,main=) h=hist(x) xfit=seq(min(x),max(x),length=100) yfit=dnorm(xfit,mean(x),sd=sd(x)) yfit=yfit*diff(h$mids[1:2])*length(x) lines(xfit, yfit, col=blue,

Re: [R] how do i plot this hist?

2010-11-09 Thread casperyc
Hi all, Thank you both. The codes work perfectly. I now use barplot(t(x.m)[-1,], names.arg = x.m[,1]) to make the 'histogram' plot. Now how do I add a dentity line to it? like 'superpose' on the graph? say if I got a Normal with mean 100 and variance 3... And I can work out the mean and

Re: [R] how do i plot this hist?

2010-11-09 Thread casperyc
Thanks! now it's just a matter of 'superposing' the density onto the barplot I am looking at the example here, but it does not seem to applicable to barplot. casper -- View this message in context: http://r.789695.n4.nabble.com/how-do-i-plot-this-hist-tp3032796p3035110.html Sent from the R

[R] how do i plot this hist?

2010-11-08 Thread casperyc
Hi all, I have the following data in abc.dat === 50 0 1 0 0 55 114 0 1 60 786 0 3 6522 324 2 3 7058 1035 1 7 7530 2568 034 80 9 293615 162 8527 216946 365

Re: [R] dot plot by group

2010-10-20 Thread casperyc
Thanks. It was 'bty=n' on a text graph, I was just having fun with it. It does not really matter that much. CasperYC -- View this message in context: http://r.789695.n4.nabble.com/dot-plot-by-group-tp2990469p3004323.html Sent from the R help mailing list archive at Nabble.com

[R] dot plot by group

2010-10-11 Thread casperyc
Hi all, I have the folloing data table %% TypeBATCH RESPONSE SHORT A 22 SHORT A 3 SHORT A 16 SHORT A 14 SHORT A 8 SHORT A 27 SHORT A 11 SHORT A 17 SHORT B 12 SHORT B 17 SHORT B 11 SHORT

Re: [R] dot plot by group

2010-10-11 Thread casperyc
Hi Spector, Yes, that is exactly what I was aiming for. Thanks. Casper -- View this message in context: http://r.789695.n4.nabble.com/dot-plot-by-group-tp2990469p2990495.html Sent from the R help mailing list archive at Nabble.com. __

Re: [R] dot plot by group

2010-10-11 Thread casperyc
And now I just wonder why the ' bty='n' ' won't work? I did dotplot(BATCH~RESPONSE,data=d,subset=Type=='SHORT',bty='n') and tried other bty parameters, none is working Casper -- View this message in context: http://r.789695.n4.nabble.com/dot-plot-by-group-tp2990469p2990500.html Sent from

[R] How to save as PDF when I used par(ask=T)

2010-04-24 Thread casperyc
== y=c(9,9,17,11,7,8,15,5) treat=c(A,C,B,C,D,A,B,D) block=c(1,1,2,2,2,3,3,3) q1.mod=aov(y~as.factor(treat)+as.factor(block)) q1.mod summary(q1.mod) par(ask=T) plot(TukeyHSD(q1.mod)) ===

[R] histogram breaks

2010-04-16 Thread casperyc
=== Q2=c( + sample(10:19,8,T), + sample(20:24,15,T), + sample(25:29,25,T), + sample(30:39,18,T), + sample(40:49,12,T), + sample(50:64,7,T), + sample(65:89,5,T) + ) hist(Q2) can give me a histogram, however, how do i get a different

[R] writing function ( 'plot' and 'if') problem

2010-04-13 Thread casperyc
=== myf=function(ds=1){ x=rnorm(10) y=rnorm(10) { #start of if if (ds==1) { list(x,y) } else (ds==2) { plot(x,y) } } # end of if } # end of function === Hi All, the problem i am having here is, that I want to be able to control the display, lf

Re: [R] writing function ( 'plot' and 'if') problem

2010-04-13 Thread casperyc
what i am triyng to do is when ds=1 give me a list ds=2 plot Thanks casper -- View this message in context: http://n4.nabble.com/writing-function-plot-and-if-problem-tp1839091p1839125.html Sent from the R help mailing list archive at Nabble.com.

[R] writing function help

2010-04-10 Thread casperyc
Hi, I have written a function, ( or 2 functions) becasue there are only tiny difference, One of them has lines == #c2 (price.mat) call option price c2=c( max(S.mat[1,3]-K,0), max(S.mat[2,3]-K,0), max(S.mat[3,3]-K,0) ) (some lines) list( Stock Value=round(S.mat,dp),

[R] writing function (plot problem)

2010-04-10 Thread casperyc
Hi, === x=rnorm(20) y=rnorm(20) t=lm(y~x) plot(t) === you will get click or hit enter to next page... how do I write a function to archieve this ? say plot(x,y) then pause, wait plot(y,x) Thanks! casper -- View this message in context:

Re: [R] R loop help

2010-04-06 Thread casperyc
Hi there, That's exactly what I want. I have checked ?combn out, but I could get the following, suppose that I want ALL possible combinations of them, as this == apply( combn(paste('x', 1:4, sep =), 2), 2, function(v) get(v[1])*get(v[2])

[R] model reparameterization

2010-04-02 Thread casperyc
== y=c(100,200,300,400,500) treatment=c(1,2,3,3,4) block=c(1,1,2,3,3) summary(lm(y~as.factor(treatment)+as.factor(block))) == The aim is to find a model that can estimate the comparison between treatment 1 with 2 and treatment 3

[R] R loop help

2010-03-26 Thread casperyc
Hi, I am tring to write a loop to compute this, == x1=c( rep(-1,4), rep(1,4) ) x2=c( rep(c(-1,-1,1,1),2) ) x3=c( rep(c(-1,1),4) ) x1*x2 x1*x3 x2*x3 suppose i have x1,x2,x3

[R] variance of discrete uniform distribution

2010-03-08 Thread casperyc
Hi all, I am REALLY confused with the variance right now. for a discrete uniform distribution on [1,12] the mean is (1+12)/2=6.5 which is ok. y=1:12 mean(y) then var(y) gives me 13 1- on http://en.wikipedia.org/wiki/Uniform_distribution_%28discrete%29 wiki the variance is

Re: [R] variance of discrete uniform distribution

2010-03-08 Thread casperyc
Hi Rolf Turner , God, it directed to the wrong page. I firstly find the formula in wiki, than tried to verify the answer in R, now, given that 143/12 ((n^2-1)/12 ) is the correct answer for a discrete uniform random variable, I am still not sure what R is calculating there? why it gives me 13?

Re: [R] test the goodness of it for negative binomial type 2

2010-03-08 Thread casperyc
Hi Achim Zeileis-4, That's very helpful. Thanks! -- View this message in context: http://n4.nabble.com/test-the-goodness-of-it-for-negative-binomial-type-2-tp1575892p1585357.html Sent from the R help mailing list archive at Nabble.com. __

[R] barplot with factors problem

2010-03-07 Thread casperyc
http://www.harding.edu/fmccown/R/#autosdatafile http://www.harding.edu/fmccown/R/#autosdatafile I am tring to get a barchat by factors, following the example in that link above. === x=c(145,40,40,120,180, 140,155,90,160,95, 195,150,205,110,160,

Re: [R] barplot with factors problem

2010-03-07 Thread casperyc
http://n4.nabble.com/file/n1583733/100307070476876317b486a941.jpg I want to get a histogram by factors. Thanks. -- View this message in context: http://n4.nabble.com/barplot-with-factors-problem-tp1583671p1583733.html Sent from the R help mailing list archive at Nabble.com.

[R] TukeyHSD model thing

2010-03-06 Thread casperyc
Hi, I am trying to reproduce a tukey test in R == x=c(145,40,40,120,180, 140,155,90,160,95, 195,150,205,110,160, 45,40,195,65,145, 195,230,115,235,225, 120,55,50,80,45 ) y2=c( rep(as.character(1),5),

[R] Fisher's LSD and tukey output thing

2010-03-06 Thread casperyc
-and how can I do that table in the attachement (page 8/9)? ( i dont know how to describe that thing, orginally produced in SAS) Q3-in the tukey output, can i ask to to show whether it's significant or not by indicating ***? (attachment page 6) Thanks! casperyc SAS output http://n4.nabble.com/file

[R] test the goodness of it for negative binomial type 2

2010-03-02 Thread casperyc
[code]library(MASS) x=c(rep(0,8096), rep(1,1629), rep(2,233), rep(3,38), rep(4,4) ) x.bar=round(mean(x),4) x.var=round(var(x),4) p.hat=round(x.bar/x.var,4) alpha.hat=round(x.bar*p.hat/(1-p.hat),4) fitdistr(x, Negative Binomial) fitdistr(x,

[R] how do I get the legend?

2010-02-20 Thread casperyc
http://en.wikipedia.org/wiki/Binomial_distribution Bino Hi there, How can I get the legend in the probability density plot in that http://en.wikipedia.org/wiki/File:Binomial_distribution_pmf.svg wiki page ? and the cumalative plot as well? # Binomial

Re: [R] legend help

2010-02-05 Thread casperyc
This is very weird!! I know the 'lty=1' command and I tried yesterday , many times!! didnt work. However, today, it worked!!! Maybe the university's computer is stupid! Thanks! -- View this message in context: http://n4.nabble.com/legend-help-tp1461992p1470216.html Sent from the R help

Re: [R] legend help

2010-02-04 Thread casperyc
Yes, that is pretty much what I want. However, there was slightly a mistake. we need to use ''rate=rate[i] and shape=shape[i] because the default is == dgamma(x, shape, rate = 1, scale = 1/rate, log = FALSE)

[R] legend help

2010-02-03 Thread casperyc
i=1 for(rate in c(2,4) ){ for(shape in c(1,3,5) ){ curve(dgamma(x,rate,shape),xlim=c(0,3),ylab=,col=i,lty=i,add=T) i=i+1 } } How can I add some legend to represent these lines? i.e. the legend is displayed as col=1 lty=1 lambda=2 theta=1 col=2 lty=2 lambda=2

Re: [R] confint for glm (general linear model)

2009-12-12 Thread casperyc
for an example, counts - c(18,17,15,20,10,20,25,13,12) outcome - gl(3,1,9); treatment - gl(3,3) glm.D93 - glm(counts ~ outcome + treatment, family=poisson()) confint(glm.D93) confint.default(glm.D93) # based on asymptotic normality to verify the confidence interval (confint.default(glm.D93))

Re: [R] confint for glm (general linear model)

2009-12-09 Thread casperyc
Zeviani - Brazil casperyc wrote: Hi, I have a glm gives summary as follows, Estimate Std. Errorz valuePr(|z|) (Intercept) -2.03693352 1.449574526 -1.405194 0.159963578 A0.01093048 0.006446256 1.695633 0.089955471 N

[R] bootstrap help

2009-12-06 Thread casperyc
Hi, I have 49 pairs in my data.frame 'data' x y 76 80 138 143 67 67 29 50 381 464 23 48 37 63 120 115 ... ... how do I get a bootstrap sample of size n=50? i have tried sample(data,size=50,replace=TRUE) and sample(data,replace=TRUE) both didnt

Re: [R] bootstrap help

2009-12-06 Thread casperyc
by the bar{x}/bar{y} ) i mis undertood that too. but now, i think what i actually need is 50 sets of this kind (49 rows) of tables from bootstrapping so that i can have 50 different 'bar{x}/bar{y}' s then to estimate the mean... Thanks. casper Ben Bolker wrote: casperyc casperyc

[R] confint for glm (general linear model)

2009-12-06 Thread casperyc
Hi, I have a glm gives summary as follows, Estimate Std. Errorz valuePr(|z|) (Intercept) -2.03693352 1.449574526 -1.405194 0.159963578 A0.01093048 0.006446256 1.695633 0.089955471 N0.41060119 0.224860819 1.826024

Re: [R] bootstrap help

2009-12-06 Thread casperyc
@r-project.org 2009/12/7 casperyc caspe...@hotmail.co.uk: Hi there, i think that's not what i was aiming for... i was aked to Generate 50 Bootstrap samples and corresponding estimates if i do data[sample(nrow(data),size=50,replace=TRUE),] it will give me a table of 50 rows ( 50 sets

Re: [R] curvedarrow (some graphics problem)

2009-06-25 Thread casperyc
Hi there, This is now the code % library(grid) vp - viewport( x = unit(0, npc), y = unit(0, npc), just = c(left, bottom), xscale = c(-1, 1) , yscale = c(-1, 1)) vp2 - viewport( # probably not needed but I had trouble placing the xaxis x = unit(0,npc), y = unit(0.5, npc), just =

Re: [R] curvedarrow (some graphics problem)

2009-06-24 Thread casperyc
baptiste auguie-4 wrote: Hi, Grid offers several functions to help drawing such graphs, see Paul Murrell's Can R Draw Graphs? (useR! 2006) I came up with this, as a quick example, vp - viewport( x = unit(0, npc), y = unit(0, npc), just = c(left, bottom), xscale = c(-1, 1)

Re: [R] curvedarrow (some graphics problem)

2009-06-24 Thread casperyc
[quote] grid.text(t, x=0.3, y=unit(-1, line), vp=vp2) Error in valid.units(units) : Invalid unit grid.text(s, 0.8, unit(-1, line), gp=gpar(col=red), vp=vp2) Error in valid.units(units) : Invalid unit[/quote] Hi, I have loaded the [grid] package. It seems that I have to load another package?

[R] curvedarrow (some graphics problem)

2009-06-22 Thread casperyc
Hi there, I just wonder how to draw this kind of picture... http://www.nabble.com/file/p24158796/b.jpg http://www.nabble.com/file/p24158796/a.jpg and this is what i have done % library(shape) library(diagram) curve(sin(x),bty=n,-8,8,yaxt=n,ylab=,xaxt=n,type=n,xlab=)