Hi Edward,

One hint to improve your code is simply stylistic. Everyone tends to

have their own favorite way, but adding spaces and line breaks can

help make code much easier to read. Also, you seemed to have used

brackets "[" instead of parentheses "(" for the first two calls to

sum(). As Berend suggested, some more clarification on what is

supposed to be happening will help you get more feedback. Lastly,

although your data is readable as is, a very convenient way to provide

the list with data is using the function, dput(). Here is a little

example demonstrating its use:

> x <- data.frame(x1 = c("c", "b", "a"),

+ x2 = factor(c("c", "b", "a")), stringsAsFactors = FALSE)

> x ## print to console, two columns appear identical

x1 x2

1 c c

2 b b

3 a a

> ## but looking at the str()ucture, they are different

> str(x)

'data.frame': 3 obs. of 2 variables:

$ x1: chr "c" "b" "a"

$ x2: Factor w/ 3 levels "a","b","c": 3 2 1

> ## use dput() to provide copy-and-pastable output

> ## that others can use to access the same data as you easy peasy

> dput(x)

structure(list(x1 = c("c", "b", "a"), x2 = structure(c(3L, 2L,

1L), .Label = c("a", "b", "c"), class = "factor")), .Names = c("x1",

"x2"), row.names = c(NA, -3L), class = "data.frame")

>

> ## now someone else is ready to go with your data

> newx <- structure(list(x1 = c("c", "b", "a"), x2 = structure(c(3L, 2L,

+ 1L), .Label = c("a", "b", "c"), class = "factor")), .Names = c("x1",

+ "x2"), row.names = c(NA, -3L), class = "data.frame")

>

> str(newx)

'data.frame': 3 obs. of 2 variables:

$ x1: chr "c" "b" "a"

$ x2: Factor w/ 3 levels "a","b","c": 3 2 1

>

And here is an example of how I might have written your function

(stylistically speaking). I am not advocating this as the "best",

"correct", or "ideal", way. Just giving another example that (I at

least) find easier to read and make sense of.

llik <- function(R_j, R_m) {

if (R_j < 0) {

sum(log(1 / (2 * pi * (sigma_j^2))) -

(1 / (2 * (sigma_j^2))*(R_j + al_j - b_j * R_m))^2)

} else if (R_j > 0) {

sum(log(1 / (2 * pi * (sigma_j^2))) -

(1 / (2 * (sigma_j^2)) * (R_j + au_j - b_j * R_m))^2)

} else if (R_j == 0) {

sum(log(pnorm(au_j, mean = b_j * R_m, sd = sigma_j) -

pnorm(al_j, mean = b_j * R_m, sd = sigma_j)))

}

}

Cheers,

Josh

On Sat, Jul 2, 2011 at 5:52 PM, EdBo <

[hidden email]> wrote:

> Hi

>

> I have developed the code below. I am worried that the parameters I want to

> be estimated are "not being found" when I ran my code. Is there a way I can

> code them so that R recognize that they should be estimated.

>

> This is the error I am getting.

>

>> out1=optim(llik,par=start.par)

> Error in pnorm(au_j, mean = b_j * R_m, sd = sigma_j) :

> object 'au_j' not found

>

> #Yet al_j,au_j,sigma_j and b_j are just estimates that balance the

> likelihood function?

>

> llik=function(R_j,R_m)

> if(R_j< 0)

> {

> sum[log(1/(2*pi*(sigma_j^2)))-(1/(2*(sigma_j^2))*(R_j+al_j-b_j*R_m))^2]

> }else if(R_j>0)

> {

> sum[log(1/(2*pi*(sigma_j^2)))-(1/(2*(sigma_j^2))*(R_j+au_j-b_j*R_m))^2]

> }else if(R_j==0)

> {

> sum(log(pnorm(au_j,mean=b_j*R_m,sd=sigma_j)-pnorm(al_j,mean=b_j*R_m,sd=sigma_j)))

> }

> start.par=c(al_j=0,au_j=0,sigma_j=0.01,b_j=1)

> out1=optim(par=start.par,llik)

>

> My Data

>

> R_j R_m

> 2e-03 0.026567295

> 3e-03 0.009798475

> 5e-02 0.008497274

> -1e-02 0.012464578

> -9e-04 0.002896023

> 9e-02 0.000879473

> 1e-02 0.003194435

> 6e-04 0.010281122

>

> Thank you in advance.

>

> Edward

> UCT

>

> --

> View this message in context:

http://r.789695.n4.nabble.com/Hint-improve-my-code-tp3641354p3641354.html> Sent from the R help mailing list archive at Nabble.com.

>

> ______________________________________________

>

[hidden email] 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.

>

--

Joshua Wiley

Ph.D. Student, Health Psychology

University of California, Los Angeles

http://www.joshuawiley.com/______________________________________________

[hidden email] mailing list

https://stat.ethz.ch/mailman/listinfo/r-helpPLEASE do read the posting guide

http://www.R-project.org/posting-guide.htmland provide commented, minimal, self-contained, reproducible code.