search for: sigmasq

...First, the optimizer returns error[8] false convergence if I estimate the functions below. I have estimated the model at first with nlm without any problems, but then I needed to add some constraints so i choose nlminb. my.loglike.normal<-function(theta){ x.start<- 1/n * sum(returns) sigmasq.start<- 1/(n-1) * sum((returns-x.start)^2) data<-c(x.start,returns) my.sigmasq<- rep(0,n+1) my.sigmasq[1]<-sigmasq.start for (i in 2:(n+1)) { my.sigmasq[i]<-(theta[1] + theta[2]*data[i-1]^2 + theta[3]*my.sigmasq[i-1])*(data[i-1]<=d1)*(my.sigmasq[i-1]&...

...d elements of returned matrix [2,3] and [3,2] are not zero as I would expect. Clearly, there is an error... What am I doing wrong? Thanks. --Dale Warning message: In if (h <= phi) { : the condition has length > 1 and only the first element will be used # function Cov.f <- function(h, sigmasq, phi) { if (h <= phi) {Cij <- sigmasq * (1 - ( 1.5 * (h/phi) - 0.5 * (h/phi)^3)) } else if (h > phi) {Cij <- 0} return(Cij) } x.coord <- c(5.7, 6.7, 9.8) y.coord <- c(42.7, 10.2, 77.4) coords <- cbind(x.coord, y.coord) distance.matrix <- as.matrix(dist(coords...

...orm(100, mean=1), rnorm(100, mean=5)) # I also tried: X <- as.matrix(c(rnorm(100, mean=1), rnorm(100, mean=5))) Xmap <- cbind(c(rep(1, 100), rep(0, 100)), c(rep(0, 100), rep(1, 100))) Xm <- mstep(modelName="V", data=X, z=Xmap) # CRASH: em(modelName=Xm$modelName, data=X, mu=Xm$mu, sigmasq=Xm$sigmasq, pro=Xm$pro) For bivariate data the following example works: # bivariate example require(mclust) Y <- as.matrix(cbind(c(rnorm(100, mean=1), rnorm(100, mean=5)), c(rnorm(100, mean=1), rnorm(100, mean=5)))) Ymap <- cbind(c(rep(1, 100), rep(0, 100)), c(rep(0, 100), rep(1, 100))) Ym...

...#### data(faithful) ################## # (2) get model ################## library(mclust) MyMixtureModel<-summary(EMclust(faithful$eruptions),faithful$eruptions) ################## # (3) plot mixture model ################## attach(MyMixtureModel) mclust1Dplot(data=faithful$eruptions,z=z,mu=mu,sigmasq=sigmasq,pro=pro,ask=FALSE,type=c("density")) do.call("mclust1Dplot",c(list(data=faithful$eruptions,ask=FALSE,type=c("density")),MyMixtureModel)) ################## # (4) plot components ################## ??? Any information you might be able to shed on this would be...

...L, 6L, 3L, 4L, 1L, 1L, 4L, 1L, 7L, 5L, 9L, 9L, 9L), .Label = c("all.all", "all.ann", "all.mon", "pft.all", "pft.ann", "pft.mon", "site.all", "site.ann", "site.mon"), class = "factor"), sigmasq = c(11430.2595455547, 12118.5387166954, 12982.4722525337, 16366.3059675243, 16650.2206047512, 19730.2121989498, 19958.3416187217, 20491.4117984889, 20647.8829877428, 21389.0300281264, 21413.7674128747, 21445.7255788782, 22002.8026436862, 22042.9802472953, 22201.0461487030, 22340...

...Dl=4.54. Using this they get a power of 0.53, but the result that I'm getting is 0.05472242. I've tried this several ways in R, but I'm not able to come up with the same number. Am I doing something wrong in the calculation of the power? Here's my code: T<-5 R<-5 sigmasq<-0.1 q<-10 Vl<-(T*R)-2 Dl<-(R*(T-1)*T*(T+1)/(12*sigmasq))*(log(1+(q/100)))^2 #Dl result is still similar power.1<-1-pf(qf(.95,(T*R-2),1,ncp=0),(T*R-2),1,ncp=Dl) Thank you for any suggestions/help. I'm using R2.4.1, on windowsXP. Erik Meesters [[alternat...

...lust package on a variable "Tuberculin indurations" recorded as mm. The file has only one variable. When I run the package I get NULL value for mu and sigma. Can anybody say why? This is the program: library("mclust") mc<-Mclust(x.trab,G=1:9,warn=TRUE) mc mc$mu sqrt(mc$sigmasq) and the output I get is > library("mclust") > mc<-Mclust(x.trab,G=1:9,warn=TRUE) There were 11 warnings (use warnings() to see them) > mc best model: unequal variance with 3 components > mc$mu NULL > sqrt(mc$sigmasq) Error in sqrt(mc$sigmasq) : Non-numeric argument...

...sing the foreach statement and I cant get it to work properly so here is what i have test<-function(){ repeat { cand2[l-1]<-rinvgamma(1,phi,lambda[l-1]) q2<-dinvgamma(cand2[l-1],phi,lambda[l-1]) p2<-cand2[l-1]^-1.5*exp(-y[l]^2/(2*cand2[l-1]))*exp(-((log(cand2[l-1])-mu_t_cand[l-1])^2)/2*sigmasq) ratio<-p2/(c[l-1]*q2) if (runif(1)< ratio) {break} } f2[l-1]<-min(p2,c[l-1]*q2) }} foreach(l=2:(n-1),.combine=c,.packages='MCMCpack',.options.nws=li st(chunkSize=250)) %dopar% test() even though i've created the array earlier called cand2 when i run this code i get the error...

Hi, I'm trying to use the library mclust for gaussian mixture on a numeric vector. The function Mclust(data,G=3) is working fine but the fitting is not optimal and is using modelNames="E". When I'm trying Mclust(data,G=3,modelName="V") I have the following message: Error in if (Sumry$G > 1) ans[c(orderedNames, "z")] else ans[orderedNames] : argument is

dear list I am student M.S. statistics in department statistics . I am working in the function "nls" in the [R 2.3.1] with 246 data and want to fit the "exp" model to vectors( v and u ) but I have a problem to use it u 5.000000e-13 2.179057e+03 6.537171e+03 1.089529e+04 1.525340e+04 1.961151e+04 2.396963e+04 2.832774e+04 3.268586e+04 3.704397e+04 4.140209e+04

I keep on receiving the message below after submitting the following line using the mclust package. m2 is a 99 X 1 column vector. * em(modelName = "E", m2, mu = c(25, 50), sigmasq=10, pro = c(0.4, 0.6)) Error in as.double.default(data) : (list) object cannot be coerced to double. Why do I receive this error? Thank, Brian C. Newquist Research Statistician Constella Health Sciences 2605 Meridian Parkway, Suite 200 Durham, NC 27713 Tel: (919) 313-7588...

Whenever I submit the following command >result <- em("E", m1, mu=c(25, 50), sigmasq=10, pro=c(0.49,0.51)) in MCLUST , I experience a problem with my Windows XP. Has anyone had this type of problem? I receive the general "send report" dialogue box and my program is no longer able to run. Should I change any of my system settings? Let me know if you think of anything...

...ond that. I have been having problems using krige.conv() to get a decent kriged map using the linear model. The code I am using from my data is as follows: >modeltest=variofit(variotest, weights=?cressie?, cov.model=?linear?, ini.cov.pars=c(80,1)) The output parameters are tausq = 9.855, sigmasq = 0.0087, phi=1.0 >krig=krige.conv(data, krige=krige.control(type.krig=?ok?, obj.model=modeltest), locations=pred.grid) At this point, krig$predict values have little to no variability (1.897 +/- 0.004), where I would expect values between 0 ? 15. By running the same data, except using an...

...y real dataset will not have peaks this well separated, but I needed to find a small example.) ################## data(faithful) library(mclust) MyMixtureModel<-summary(EMclust(faithful$eruptions),faithful$eruptions) attach(MyMixtureModel) mclust1Dplot(data=faithful$eruptions,z=z,mu=mu,sigmasq=sigmasq,pro=pro,ask=FALSE,type=c("density")) do.call("mclust1Dplot",c(list(data=faithful$eruptions,ask=FALSE,type=c("density")),MyMixtureModel)) # plot components ??? ################## Any information you might be able to shed on this would be very much apprec...

.... The user defines an arbitrary set of parameters as the function's arguments, which are used to construct this covariance matrix. Within the C code these parameters and hence covariance matrix are iteratively updated. The user's function might look like this: cov.function <- function(sigmasq, tausq, phi){ sigmasq*exp(-phi*D) + tausq*diag(n) } The way my C code is now, the iteratively proposed values for the cov.fun parameters are in the same order as the parameters used in the cov.fun function (i.e., same order as the formal list). So what I'd like to do is just set the value p...

...thin = 100) #trial error outmcmc.5 = glsm.mcmc(geoData9093, model= model.5, mcmc.input = mcmc.5) mcmcobj.5 = prepare.likfit.glsm(outmcmc.5) lik.5 = likfit.glsm(mcmcobj.5, ini.phi = 0.3, fix.nugget.rel = F)/ And the summary of lik.5 is: likfit.glsm: estimated model parameters: beta sigmasq phi tausq.rel "1.2781" "0.5193" "0.0977" "0.0069" likfit.glsm : maximised log-likelihood = 43.62 I'm fairly new to geostatistics, but I thought using a factor variable as covariable would give me 4 intercepts (beta) as I have 4 levels in my...

...)).? However, when I use rnorm to generate some fake data, mclust works consistently well.? As I am quite new to R, I was wondering whether anyone could help me out.? I am using the em function in mclust, with parameters that are arbitrarily set Pro <- c(0.33,0.33,0.33) Mean<- c(-0.5,0,0.5) Sigmasq <- c(0.05,0.05,0.05) And a sample of my data set is c(-0.016133,-0.07668,-0.000625,0.031329,-0.011094,0.014199,-0.014141,0.036836,-0.007695,-0.011738,-0.021953,0.046565,0.010859,-0.050313,-0.167813,0.01543,-0.057598,-0.034336,0.182275,0.032959,0.01918,0.009248,-0.195273,-0.00918,-0.017813,0.00...

...1 - log(1+exp(eta1)); eta2 = a0 + a1*lnndsTs; llnorm = -1/(2*sigsq)*(lnbldT - eta2)**2 - .5*log(sigsq); ll = llbin + llnorm; end; else do; eta2 = a0 + a1*lnndsTs; eta1 = b0 + b1*eta2 + u; llbin = anybc.cens.ind*eta1 - log(1+exp(eta1)); ll = llbin; end; sigma2 = sigmasq*b1**2 /*variance of random effect; model anybc.cens.indic ~ general(ll); random u ~ normal(0,sigma2) subject = CaseID; */; run; #################################################################### R Code with error message: me.km.nlme <- nlme(model = anybc.cens.indic ~ vs.flag*((anybc.cens.ind...

...gramming error due to my inexperience with R. I hope some fo you can show me where the error is. Thanks Ruben mcolasim7<-function(N.sim,pcell){ # Lattice definition x<-seq(1,180,1) y<-seq(1,540,1) #Gaussian process param<-c(6.63,2.24,1.82,4.36) names(param)<-c("beta","sigmasq","tausq","varphi") #Simulations loop - Function GaussRF from package RandomFields for(i in 1:N.sim){ mcola<-GaussRF(x=x,y=y,param=param,grid=TRUE,model="gauss") #Process thininng out with 'pcell': probability that the process will manifest itself for(j...

...ve a problem in fitting lines() of the normal distributions identified with Mclust on a histogram or a mclust1Dplot. Here is some sample code to explain : set.seed(22) foo <- c(rnorm(400, 10, 2), rnorm(500, 17, 4)) mcl <- Mclust(foo, G=2) mcl.sd <- sqrt(mcl$parameters$variance$sigmasq) mcl.size <- c(length(mcl$classification[mcl$classification==2]), length(mcl$classification[mcl$classification==1])) x <- pretty(c(0:44), 100) #### my plot of histogram and lines of normal distributions #### SEEMS OK (or am I wrong ?) using frequencies : histA <- hist(foo, bre...