+ Reply to Thread
Results 1 to 5 of 5

Thread: Error in R-code for evaluation of optimal design

  1. #1
    Points: 8, Level: 1
    Level completed: 15%, Points required for next Level: 42

    Posts
    2
    Thanks
    0
    Thanked 0 Times in 0 Posts

    Error in R-code for evaluation of optimal design




    For my project I have to choose an optimal design for a given problem and I have to evaluate this design. The problem follows a linear regression model (cubic relationship), with a variance of 23. Out of the restrictions that were given, I know the 4 coefficients of the parameters.

    My evaluation (R-code below) shows biased results and there would be no correct coverage of my 95% CI. Where did I go wrong in my script?


    # de model parameters (data generend model)
    beta0<-10
    beta1<-69.8333333
    beta2<- -15.7500000
    beta3<-0.9166667

    # effect dat moet geschat worden
    effect<-beta1

    #steekproefgroottes
    n<-8 # steekproefgrootte is 8 voor een G-optimal design

    #initialisatie van vectoren
    beta1.N<-c()
    s2.beta1.N<-c()
    cnt<-0

    #dataframe voor component
    efficientie<-data.frame(component=c(10,0,2.45,3.12,7.25,10,0,7.25))

    #aantal simulaties
    N<-10000
    for(i in 1:N){
    X1<-c(10,0,2.45,3.12,7.25,10,0,7.25)
    Y<-beta0+beta1*X1+beta2*X1^2+beta3*X1^3

    #steek alle data in de data matrix
    efficientie$component<-c(10,0,2.45,3.12,7.25,10,0,7.25)
    efficientie$component2<-X1^2
    efficientie$component3<-X1^3
    efficientie$score<-Y

    #fit het werkingsmodel aan de data
    model.fit<-lm(score~component+component2+component3,data=efficientie)

    #extraheer de parameterschatting beta1 uit het gefitte model
    beta1.N[i]<-summary(model.fit)$coeff[2,1]
    #extraheer de geschatte variantie op beta1_hat
    s2.beta1.N[i]<-summary(model.fit)$coeff[2,2]^2
    #95% betrouwbaarheidsinterval
    CI<-confint(model.fit)
    #Nagaan van de 95% bedekking
    if((effect>CI[2,1])&(effect<CI[2,2]))cnt<-cnt+1
    }


    mean(beta1.N)

    var(beta1.N)
    mean(s2.beta1.N)

    cnt/N

  2. #2
    TS Contributor
    Points: 12,030, Level: 71
    Level completed: 95%, Points required for next Level: 20
    rogojel's Avatar
    Location
    I work in Europe, live in Hungary
    Posts
    1,451
    Thanks
    159
    Thanked 330 Times in 310 Posts

    Re: Error in R-code for evaluation of optimal design

    Hi,
    I may not understand the code fully, but it seems to me that you are calculating N=10000 times the exactly same model, with no random variations whatsoever. Do I miss something?

    regards

  3. #3
    Points: 8, Level: 1
    Level completed: 15%, Points required for next Level: 42

    Posts
    2
    Thanks
    0
    Thanked 0 Times in 0 Posts

    Re: Error in R-code for evaluation of optimal design

    Hi,

    Thanks for your answer! Today I also concluded that this lus is always calculating the same model. But I'm still stuck, I don't know with what I have to compare this model with. With random Y-calculations? But then I get a strong biased beta1 estimator.
    I have to do an evaluation of the chosen design (G-optimal design) and the point estimators. A cubic relationship is assumed. So what I have done now, is calculating the Y-values with the cubic relationship, the optimal designpoints, and the calculated coefficients (bčtas). But indeed, it's not smart to do estimations of beta with this model.. But at the moment, I have no idea of how to evaluate this.

  4. #4
    TS Contributor
    Points: 12,030, Level: 71
    Level completed: 95%, Points required for next Level: 20
    rogojel's Avatar
    Location
    I work in Europe, live in Hungary
    Posts
    1,451
    Thanks
    159
    Thanked 330 Times in 310 Posts

    Re: Error in R-code for evaluation of optimal design

    What is the question you want to answer with this code?

  5. #5
    Human
    Points: 12,596, Level: 73
    Level completed: 37%, Points required for next Level: 254
    Awards:
    Master Tagger
    GretaGarbo's Avatar
    Posts
    1,352
    Thanks
    455
    Thanked 461 Times in 401 Posts

    Re: Error in R-code for evaluation of optimal design


    We are not so good in understanding Dutch. Your code explanations are cryptic.
    (Some of us are not so good in English either but we understand the elementary part.)

    "Today I also concluded that this lus is always...."
    It seems like "lus" in Dutch means "loop".

    Make a fully reproducible example and explain what the problem is and highlight the code and click on the # symbol.

+ Reply to Thread

           




Posting Permissions

  • You may not post new threads
  • You may not post replies
  • You may not post attachments
  • You may not edit your posts






Advertise on Talk Stats