search for: yhat

...his problem I have by running the following example using my code at the end of the email: spp.dat <- matrix(c(4,2,0,2,0,5,1,3,2,4,0,2,2,0,3,1),4,4) env.dat <- matrix(c(1.5,2.3,2,1.6,0.9,0.8,1.2,1.5),4,2) test.rda <- rda(spp.dat, env.dat) Results in the following error: Error in df %*% t(Yhat) : non-conformable arguments This is referring to the part of the code where I wish to do the following: S_Yhat'Yhat = [1/(n-1)] Yhat' Yhat, where [1/(n-1)] = df = degrees of freedom. is effectively dividing my t(Yhat) %*% Yhat by the degrees of freedom, n being th...

I have a question about passing arguments to the function f that nlm minimizes. I have no problems if I do this: x<-seq(0,1,.1) y<-1.1*x + (1-1.1) + rnorm(length(x),0,.1) fn<-function(p) { yhat<-p*x+(1-p) sum((y-yhat)^2) } out<-nlm(fn,p=1.5,hessian=TRUE) But I would like to define fn<-function(x,y,p) { yhat<-p*x+(1-p) sum((y-yhat)^2) } so that I can pass the x and y variables along, doing something like out<-nlm(fn,x=this, y=that,p=1.5,hes...

...ute the Hosmer-Lemeshow 'goodness-of-fit' test cd.full_model = glm(formula = Collaterals ~ CHF + Age + CABG + relevel (as.factor (num.obst.vessels),"one") + Current.smoker + DM + HTN + ace.inhibitor + MI, family = binomial(link = "logit")) hosmerlem = function(y, yhat, g=10) { cutyhat = cut(yhat, breaks = quantile(yhat, probs=seq(0, 1, 1/g)), include.lowest=TRUE) obs = xtabs(cbind(1 - y, y) ~ cutyhat) expect = xtabs(cbind(1 - yhat, yhat) ~ cutyhat) chisq = sum((obs - expect)^2/expect) P = 1 - pchisq(chisq, g - 2) return(list(chisq=chisq,p...

Hi, See below I reply your message for <https://stat.ethz.ch/pipermail/r-help/2008-April/160966.html>[R] predict.glm & newdata posted on Fri Apr 4 21:02:24 CEST 2008 You say it ##works fine but it does not: if you look at the length of yhat2, you will find 100 and not 200 as expected. In fact predict(reg1, data=x2) gives the same results as predict(reg1). So I am still looking for a solution of this problem. best regards Maurice =========== Hi, could it be the "newdata" argument? When I run predict with the newdata...

...he */retain/ *statement. For people that want to read complicated explanations to help me: I have a system of two equations written as difference equations here. To boil it down. I have a dataframe with three variables y, X1, X2 which are measured. Now, I want to estimate a model which says y= yhat+epsillon as in standard nonlinear regression. Where yhat is the predicted value of y. yhat can be calculated as follows: # This code illustrates how I could simulate the expected values of y if I knew the values of the parameters tau and b. # but in reality I would like to estimate them. # code...

...is PDF: http://www.utdallas.edu/~herve/abdi-prc-pretty.pdf I created a small CSV called, students.csv that contains the following data: s1 14 4 1 s2 23 4 2 s3 30 7 2 s4 50 7 4 s5 39 10 3 s6 67 10 6 Col headers: Student id, Memory span(Y), age(X1), speech rate(X2) Now the expected results are: yHat[0]:15.166666666666668 yHat[1]:24.666666666666668 yHat[2]:27.666666666666664 yHat[3]:46.666666666666664 yHat[4]:40.166666666666664 yHat[5]:68.66666666666667 This is based on the following equation (given in the PDF): Y = 1.67 + X1 + 9.50 X2 I ran the following commands in R: data = read.table(&q...

...eak' is available at the following URL        http://www.stat.sc.edu/~hitchcock/diseaseoutbreakdata.txt The R code as mentioned therein is #### # A function to do the Hosmer-Lemeshow test in R. # R Function is due to Peter D. M. Macdonald, McMaster University. # hosmerlem <- function (y, yhat, g = 10) { cutyhat <- cut(yhat, breaks = quantile(yhat, probs = seq(0, 1, 1/g)), include.lowest = T) obs <- xtabs(cbind(1 - y, y) ~ cutyhat) expect <- xtabs(cbind(1 - yhat, yhat) ~ cutyhat) chisq <- sum((obs - expect)^2/expect) P <- 1 - pchisq(chisq, g -...

...# set -ve Var(Y) to zero #+++++++++++++++++++++++++++ Fit WLS +++++++++++++++++++++++++++++++++++++++++++ fitwls <- lm(formula = fmla1,weights = Weights,data = data.frame(cbind(X,Y))) bhat <- coef(fitwls) ############################## Y = Log(Z) Scale #################################### Yhat <- X%*%bhat # predicted values mu <- mean(Y) To <- Y - mu Er <- Y - Yhat Re <- Yhat - mu lgSST <- sum(Weights*(To)^2) # log SST lgSSE <- sum(Weights*(Er)^2) # log SSE lgSSR <- sum(Weights*(Re)^2) # log SSR lgR-sq &...

I am trying to fit a nonlinear regression model to data. There are several predictor variables and 8 parameters. I will write the model as Y ~ Yhat(theta1,...,theta8) OK, I can do this using nls() - but "only just" as there are not as many observations as might be desired. Now the problem is that we have a factor "Site" and I want to include a corresponding error component. I tried something like (excuse the doctored R...

...check to see if we fed in the right kind of object stopifnot(family(obj)$family=="binomial" && family(obj)$link=="logit") y = obj$model[[1]] # the double bracket (above) gets the index of items within an object if (is.factor(y)) y = as.numeric(y)==2 yhat = obj$fitted.values cutyhat=cut(yhat,quantile(yhat,0:g/g),include.lowest=TRUE) obs = xtabs(cbind(1 - y, y) ~ cutyhat) expect = xtabs(cbind(1 - yhat, yhat) ~ cutyhat) if (any(expect < 5)) # warning("Some expected counts are less than 5. Use smaller number of groups")...

...is calculated based on the fixed and pre-determined cut-off points of the estimated probability of success. One of them is calculated based on the percentiles of estimated probabilities. In the previous post,i find that the Hosmer and Lemeshow?s test how to use in R. hosmerlem <- function (y, yhat, g = 10) { cutyhat <- cut(x, breaks = quantile(yhat, probs = seq(0, 1, 1/g)), include.lowest = T) obs <- xtabs(cbind(1 - y, y) ~ cutyhat) expect <- xtabs(cbind(1 - yhat, yhat) ~ cutyhat) chisq <- sum((obs - expect)^2/expect) P <- 1 - pchisq(chisq, g - 2)...

Hello, I have the following function that receives a "function pointer" formal parameter name "fnc": loocv <- function(data, fnc) { n <- length(data.x) score <- 0 for (i in 1:n) { x_i <- data.x[-i] y_i <- data.y[-i] yhat <- fnc(x=x_i,y=y_i) score <- score + (y_i - yhat)^2 } score <- score/n return(score) } I would like to use it like this: ## ## Estimator function using Gaussian Kernel ## gaussiankernel <- function(x,y,h) { modelks <- ksmooth(x,y,kernel=&qu...

Hi all - I'm stumped by the following mdl <- glm(resp ~ . , data = df, family=binomial, offset = ofst) WORKS yhat <- predict(mdl) WORKS yhat <- predict(mdl,newdata = df) FAILS Error in drop(X[, piv, drop = FALSE] %*% beta[piv]) : subscript out of bounds I've tried without offset, quoting binomial. The offset variable ofst IS in df. Previous postings indicate possible names(df) problems (renami...

...1000) e1<-rnorm(1000)/3 e2<-rnorm(1000)/3 e3<-rnorm(1000)/3 y<-1-(1-pnorm(-2+0.33*x1+0.66*x2+1*x3+e1)*1-(pnorm(1+1.5*x4-0.25*x5+e2)*pnorm(1+0.2*x6+e3))) y <- y>runif(1000) dat<-data.frame(y = y, x1 = x1, x2 = x2, x3 = x3) g<-glm(y~., data = dat, family = binomial) summary(g) yhat<-predict(g, dat) Call: glm(formula = y ~ ., family = binomial, data = dat) Deviance Residuals: Min 1Q Median 3Q Max -1.8383 -1.3519 0.7638 0.9249 1.3698 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) 0.71749 0.06901 1...

...it = 50, trace = F))   > # Option 1: Hosmer-Lemeshow test > mod.fit <- glm(formula = no.NA$repcnd ~  no.NA$svl, family = binomial(link = logit), data =  no.NA, na.action = na.exclude, control = list(epsilon = 0.0001, maxit = 50, trace = F))                >  hosmerlem <- function (y, yhat, g = 10) { cutyhat <- cut(yhat, breaks = quantile(yhat, probs = seq(0, 1, 1/g)), include.lowest = T) obs <- xtabs(cbind(1 - y, y) ~ cutyhat) expect <- xtabs(cbind(1 - yhat, yhat) ~ cutyhat) chisq <- sum((obs - expect)^2/expect) P <- 1 - pchisq(chisq, g - 2) c("X^2" = chisq...

...Residual deviance: 285.04  on 216  degrees of freedom   (8 observations deleted due to missingness) AIC: 289.04   Number of Fisher Scoring iterations: 3 ------------------------------------------------------------------------------- > Hosmer-Lemeshow test > > hosmerlem <- function (y, yhat, g = 10) + { + cutyhat <- cut(yhat, breaks = quantile(yhat, probs = seq(0, 1, 1/g)), include.lowest = T) + obs <- xtabs(cbind(1 - y, y) ~ cutyhat) + expect <- xtabs(cbind(1 - yhat, yhat) ~ cutyhat) +  chisq <- sum((obs - expect)^2/expect) + P <- 1 - pchisq(chisq, g - 2) + c("X^...

....seed(121) x1=.04 for (i in 1:14) x1[i+1]=x1[i]*(1+rnorm(1)*.008)+.00025 x2=.08 for (i in 1:14) x2[i+1]=x2[i]*(1+rnorm(1)*.03)-.0018 x3=.01 for (i in 1:14) x3[i+1]=x3[i]*(1+rnorm(1)*.15)-.0008 b=matrix(c(0.6,0.0,0.4)) x=matrix(cbind(x1,x2,x3),ncol=3) y=x%*%b # the 'real' y yhat=y+runif(15)*.006 # the observed y plot(x=1:15,ylim=c(min(x1,x2,x3),max(x1,x2,x3))) matlines(cbind(x,y,yhat)) # Add a constant to x (for alpha) x=cbind(1,x) # "normalization" fun to make the rest of the x's add up to 1 normalize=function(x)c(x[1], x[2:length(x)]/sum(x[2:length(x)]...

Dear All I am trying to estimate a mixture of regression and get the following error using the mixreg package: Error in y - yhat : non-conformable arrays The instruction I used were: x <- as.matrix(LRHUN) y <- as.matrix(LRINTER) TS <- list(list(beta=c(3.0,1.0),sigsq=1,lambda=0.4), list(beta=c(0.0,1.0),sigsq=1,lambda=0.6)) prova <- mixreg(x,y, ncomp=2, theta.start=TS) The data set has about 1200 o...

...e this: # need this library(rms) # setup data d <- data.frame(x=rnorm(100), y=rnorm(100)) dd <- datadist(d) options(datadist='dd') # fit model l <- ols(y ~ rcs(x), data=d) # predict along original limits of data l.pred <- Predict(l) # plot of fit and conf. int. xYplot(Cbind(yhat, lower, upper) ~ x, data=l.pred, method='filled', col.fill=grey(0.9), col=1, type='l') Is there any way in which I can turn this figure on its side, by plotting x ~ y... something like this: # doesn't work xYplot(x ~ Cbind(yhat, lower, upper), data=l.pred, method='filled...

...to use the GLS function to fit the linear regression with the correlation structure AR(3) for the response (or residual). I hope the residuals after the GLS regression will be independent judged by Box-Ljung test. However, I dont know how the residuals are defined in the GLS function. Is it just y-yhat or y-yhat times the (covariance matrix)^(-1/2). Because y-yhat still has the AR(3) covariance structure and surely be rejected by the Box-Ljung test. The latter will be independent if the assumption of AR(3) correlation structure is right. Any suggestion are highly appreciated. Thx!