Hi all,
I am facing difficulty on how to use bootstrap sampling and
below is my example of function.
Read a data , use some functions and use iteration to find the solution(
ie, convergence is reached). I want to use bootstrap approach to do it
several times (200 or 300 times) this whole process and see the
distribution of parameter of interest.
Below is a small example that resembles my problem. However, I found out
all samples are the same. So I would appreciate your help on this case.
#**************************************
rm(list=ls())
xx <- read.table(textConnection(" y x
11 5.16
11 4.04
14 3.85
19 5.68
4 1.26
23 7.89
15 4.25
17 3.94
7 2.35
17 4.74
14 5.49
11 4.12
17 5.92"), header=TRUE)
data <- as.matrix(xx)
closeAllconnections()
Nt <- NULL
for (Ncount in 1:100)
{
y <- data[,1]
x <- data[,2]
n <- length(x)
X <- cbind(rep(1,n),x) #covariate/design matrix
obeta<- c(1,1) #previous/starting values of beta
nbeta <- c(0,0) #new beta
iter=0
while(crossprod(obeta-nbeta)>10^(-12))
{
nbeta <- obeta
eta <- X%*%nbeta
mu <- eta
mu1 <- 1/eta
W <- diag(as.vector(mu1))
Z <- X%*%nbeta+(y-mu)
XWX <- t(X)%*%W%*%X
XWZ <- t(X)%*%W%*%Z
Cov <- solve(XWX)
obeta <- Cov%*%XWZ
iter <- iter+1
cat("Iteration # and beta1= ",iter, nbeta, "\n")
}
Nt[Ncount] <- nbeta[1,1]
}
Nt
summary(Nt)
#**************e*****************************************
[[alternative HTML version deleted]]
I didn't see bootstrap steps in your code. One way to modify your codes
for (Ncount in 1:100)
{
b.data<-data[sample(1:nrow(data),replace=T),]
y <-b.data[,1]
x <-b.data[,2]
n <- length(x)
... ### make appropriate changes if needed
}
Weidong Gu
On Wed, Jul 20, 2011 at 6:09 PM, Val <valkremk at gmail.com>
wrote:> Hi all,
>
> I am facing difficulty on ?how to use bootstrap sampling and
> below is my example of function.
>
> Read a data , use some functions and ?use iteration to find the solution(
> ie, convergence is reached). ?I want to use bootstrap approach to do it
> several times (200 or 300 times) this whole process ?and see the
> distribution of parameter of interest.
>
> Below is a small example that resembles my problem. However, ?I ?found out
> all samples are the same. So I would appreciate your help on this case.
>
> #**************************************
> rm(list=ls())
> ?xx <- read.table(textConnection(" y x
> ? ?11 5.16
> ? ?11 4.04
> ? ?14 3.85
> ? ?19 5.68
> ? ?4 1.26
> ? ?23 ?7.89
> ? ?15 4.25
> ? ?17 3.94
> ? ?7 2.35
> ? ?17 4.74
> ? ?14 5.49
> ? ?11 4.12
> ? ?17 5.92"), header=TRUE)
> ? ?data <- as.matrix(xx)
> ? ?closeAllconnections()
>
> Nt <- NULL
> for (Ncount in 1:100)
> ?{
> ? ?y <- data[,1]
> ? ?x <- data[,2]
> ? ?n <- length(x)
>
> ? ?X <- cbind(rep(1,n),x) ? ? ? ? ? ? ? ? #covariate/design matrix
> ? ?obeta<- c(1,1) ? ? ? ? ? ? ? ? ? ? ? ? #previous/starting values of
beta
>
> ? ?nbeta <- c(0,0) ? ? ? ? ? ? ? ? ? ? ? ?#new beta
> ? ?iter=0
>
> ?while(crossprod(obeta-nbeta)>10^(-12))
> ? {
> ? ?nbeta <- obeta
> ? ?eta ? <- X%*%nbeta
> ? ?mu ? ?<- eta
> ? ?mu1 ? <- 1/eta
> ? ?W ? ? <- diag(as.vector(mu1))
> ? ?Z ? ? <- X%*%nbeta+(y-mu)
> ? ?XWX ? <- t(X)%*%W%*%X
> ? ?XWZ ? <- t(X)%*%W%*%Z
> ? ?Cov ? <- solve(XWX)
> ? ?obeta <- Cov%*%XWZ
> ? ?iter ?<- iter+1
>
> ? ?cat("Iteration # ?and beta1= ",iter, nbeta, "\n")
> ? ?}
>
> ?Nt[Ncount] <- nbeta[1,1]
> }
> Nt
> summary(Nt)
> #**************e*****************************************
>
> ? ? ? ?[[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org 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.
>
In order to apply the bootstrap, you must resample, uniformly at random
from the independent units of measurement in your data. Assuming that
these represent the rows of 'data', consider the following:
est <- function(y, x, obeta = c(1,1), verbose=FALSE) {
n <- length(x)
X <- cbind(rep(1, n), x)
nbeta <- c(0,0)
iter <- 0
while(crossprod(obeta-nbeta)>10^(-12)) {
nbeta <- obeta
eta <- X%*%nbeta
mu <- eta
mu1 <- 1/eta
W <- diag(as.vector(mu1))
Z <- X%*%nbeta+(y-mu)
XWX <- t(X)%*%W%*%X
XWZ <- t(X)%*%W%*%Z
Cov <- solve(XWX)
obeta <- Cov%*%XWZ
iter <- iter+1
if(verbose)
cat("Iteration # and beta1= ",iter, nbeta,
"\n")
}
return(nbeta[1,1])
}
boot <- function(data, reps) {
n <- nrow(data)
Nt <- vector('numeric', length=reps)
for(Ncount in 1:reps) {
#resample the rows of data
bdata <- data[sample(1:n,n,replace=TRUE),]
#recompute and store estimate
Nt[Ncount] <- est(bdata[,1], bdata[,2])
}
return(Nt)
}
stem(boot(data,1000),width=60)
The decimal point is at the |
-3 | 4
-2 |
-1 | 2
-0 | 888666666555444333333333332222222222211111111111
0 | 000000000000001111111111111111111111111111222222+400
1 | 000000000000000000000000000000011111111111111111+203
2 | 000000000000000011111111222222222223333444444444+23
3 | 000000001111111111222222222223344444445555556666
4 | 11222233445556666789
5 | 02334446677899
6 | 1112334455778
7 | 000011235568
8 | 001799
9 | 0259
10 | 1446
11 | 19
12 | 48
13 | 8
14 | 024
15 |
16 |
17 | 0788
18 |
19 | 1
On Wed, 2011-07-20 at 18:09 -0400, Val wrote:> Hi all,
>
> I am facing difficulty on how to use bootstrap sampling and
> below is my example of function.
>
> Read a data , use some functions and use iteration to find the solution(
> ie, convergence is reached). I want to use bootstrap approach to do it
> several times (200 or 300 times) this whole process and see the
> distribution of parameter of interest.
>
> Below is a small example that resembles my problem. However, I found out
> all samples are the same. So I would appreciate your help on this case.
>
> #**************************************
> rm(list=ls())
> xx <- read.table(textConnection(" y x
> 11 5.16
> 11 4.04
> 14 3.85
> 19 5.68
> 4 1.26
> 23 7.89
> 15 4.25
> 17 3.94
> 7 2.35
> 17 4.74
> 14 5.49
> 11 4.12
> 17 5.92"), header=TRUE)
> data <- as.matrix(xx)
> closeAllconnections()
>
> Nt <- NULL
> for (Ncount in 1:100)
> {
> y <- data[,1]
> x <- data[,2]
> n <- length(x)
>
> X <- cbind(rep(1,n),x) #covariate/design matrix
> obeta<- c(1,1) #previous/starting values of
beta
>
> nbeta <- c(0,0) #new beta
> iter=0
>
> while(crossprod(obeta-nbeta)>10^(-12))
> {
> nbeta <- obeta
> eta <- X%*%nbeta
> mu <- eta
> mu1 <- 1/eta
> W <- diag(as.vector(mu1))
> Z <- X%*%nbeta+(y-mu)
> XWX <- t(X)%*%W%*%X
> XWZ <- t(X)%*%W%*%Z
> Cov <- solve(XWX)
> obeta <- Cov%*%XWZ
> iter <- iter+1
>
> cat("Iteration # and beta1= ",iter, nbeta, "\n")
> }
>
> Nt[Ncount] <- nbeta[1,1]
> }
> Nt
> summary(Nt)
> #**************e*****************************************
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org 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.