similar to: Vectorize (scalar) function

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]) k=k+1 } } Basically, I want to investigate for what values of mu and sigma, the integral is divergent. Is there another way

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 stopped working. On the other hand, Maple gives me a value of 0.5005299403. It is an

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 ----- ################################################### PhD candidate in Statistics School

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, it should be able to use 'predict' and 'plot' and other

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", lwd=2) } ############################# individually, it worked fine however, if I used par(mfrow=c(2,2))

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 http://stat.ethz.ch/R-manual/R-devel/library/utils/html/prompt.html created the file plot.default.rd in my directory. But

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: http://r.789695.n4.nabble.com/how-to-add-frequencies-to-barplot-tp3051923p3051923.html Sent from the R help

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 =1) Thanks. Casper -- View this message in context:

=========================== 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 ds=1, i want to just have a list, but it seem to always plot... Thanks. casper -- View this

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, it's not working... Every time, it goes to the end, the only thing I can click is

Hi, I have a glm gives summary as follows, Estimate Std. Error z value Pr(>|z|) (Intercept) -2.03693352 1.449574526 -1.405194 0.159963578 A 0.01093048 0.006446256 1.695633 0.089955471 N 0.41060119 0.224860819 1.826024 0.067846690 S -0.20651005 0.067698863 -3.050421 0.002285206 then I use confint(k.glm)

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 dat3.res dat4.res dat5.res stored in the workspace. This is only an example of a complex

Hi all, I have the following data in abc.dat ======================= 50 0 1 0 0 55 1 14 0 1 60 7 86 0 3 65 22 324 2 3 70 58 1035 1 7 75 30 2568 0 34 80 9 2936 15 162 85 27 2169 46 365 90 80 1439 212 432 95 236 1670 521 281 100 332 827 709 172 105 156 311 556

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, 45,40,195,65,145, 195,230,115,235,225, 120,55,50,80,45 ) y2=c( rep(as.character(1),5), rep(as.character(2),5),

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 values, say one of the following: k=20 mu=-1.978295 casigma=0.008326927 >

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 like to achieve more. Apart from using a loop, is there are fast way to 'add' the titles to be

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: http://n4.nabble.com/writing-function-plot-problem-tp1835723p1835723.html Sent from the R help mailing

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', xlab="",ylab="" ) lines(p,g1,lty=1,col=1) lines(p,g2,lty=1,col=2)

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 theta=3 col=3 lty=3 lambda=2 theta=5 col=4 lty=4 lambda=4 theta=1 col=5 lty=5 lambda=4 theta=3 col=6 lty=6 lambda=4 theta=5

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 (j<4){ if (j==2) {x=x+c(-1,1)*0.5} ci[j,]= x+qnorm(1-alpha/2)^2/2+ c(-1,1)*qnorm(1-alpha/2)* sqrt(x+qnorm(1-alpha/2)^2/4) j=j+1 if (j==3) { # exact x=x-c(-1,1)*0.5