[R] Optimize function in R: unable to find maximum of logistic function

2013-10-28 Thread Shantanu MULLICK
Hello Everyone,

I want to perform a 1-D optimization by using the optimize() function. I
want to find the maximum value of a logistic function. The optimize()
function gives the wrong result.

My code:
f= function (k) {
T_s = 20
result = (2- 2/(1+ exp(-2*T_s*k)))
return(result)
}
optimize(f, c(0, 5), tol = 0.1, maximum= TRUE)

The maximum value for the function happens at k=0, and the maximum value is
1. Yet  the optimise function, says that the maximum value happens at k=
4.9995, and the maximum value is 0.

Thanks in advance!

Warm Regards,
Shantanu

[[alternative HTML version deleted]]

__
R-help@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.


Re: [R] Optimize function in R: unable to find maximum of logistic function

2013-10-28 Thread Shantanu MULLICK
Hello,

Thanks a lot for the replies.

I tried to use the optim function in R, and chose the L-BFGS-B method
of optimization. Now the result is very sensitive to the starting values. I
use the same function, but I give the line of  code I use for optim.

optim(0.25, f, method= L-BFGS-B, lower = 0, upper=5). Here 0.25 is the
starting value.

Result: $ par = 0, $Value = 1 (this is correct). $par indicates that the
maximum value of the function is at 0, and $ Value says that the maximum
value is 1.

However, I produced a quick table for the results with different starting
values:
Starting Value $par $Value
 Result
0.25 0  1
  Correct
0.5  0.4999  4.1223e-09
Wrong
11   0
   Wrong
3 3   0
  Wrong

I wanted to optimise the logistic function, as the logistic function
approximates an indicator function. The indicator function is defined as
follows:  I(u), where I(u) = 1, where u = 0, and I(u) = 0, for all other
 u. u is strictly positive.

My objective function contains a indicator function.As the indicator
function is non-continuous, the optimization would become non-convex. To
bypass this limitation, I used the logistic function to approximate the
indicator function.

A simple example: f(u) = I(u) + u/6;Maximise f(u) where u is between (0
to 5).
Correct Answer: f(u) =1, where u=0.
Code:
f= function (k) {
T_s = 20
result = (2- 2/(1+ exp(-2*T_s*k))) + (k/6)
return(result)
}
optim(0.15, f, method= L-BFGS-B, lower = 0, upper=5,control =
list(fnscale=-1) )
## Here 0.15 is the starting value. I produce a similar table as above:

Starting Value $par $Value
 Result
0.15 0  1
  Correct
0.16 5   0.833
 Wrong
0.17 5   0.833
 Wrong
2   5   0.833
   Wrong

Questions:
1. If my objective function contains an Indicator function, is there an
efficient way to implement it?
2. Does R have any quick non-convex optimiser functions?

Best Regards,
Shantanu Mullick
PhD Candidate
ESSEC Business School

On 28 October 2013 22:01, Rolf Turner r.tur...@auckland.ac.nz wrote:


 This could be described as a bug, perhaps.  Or it could be described as
 an indication that numerical optimization is inevitably tricky. Notice that
 if you narrow down your search interval from [0,5] to [0,0.5] you get the
 right answer:

  optimize(f, c(0, 0.5), maximum= TRUE,tol=1e-10)
 $maximum
 [1] 4.192436e-11

 $objective
 [1] 1

 I guess there's a problem with finding a gradient that is effectively
 (numerically) zero when k is equal to 5.


 cheers,

 Rolf Turner



 On 10/29/13 06:00, Shantanu MULLICK wrote:

 Hello Everyone,

 I want to perform a 1-D optimization by using the optimize() function. I
 want to find the maximum value of a logistic function. The optimize()
 function gives the wrong result.

 My code:
 f= function (k) {
 T_s = 20
 result = (2- 2/(1+ exp(-2*T_s*k)))
 return(result)
 }
 optimize(f, c(0, 5), tol = 0.1, maximum= TRUE)

 The maximum value for the function happens at k=0, and the maximum value
 is
 1. Yet  the optimise function, says that the maximum value happens at k=
 4.9995, and the maximum value is 0.



[[alternative HTML version deleted]]

__
R-help@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.


[R] R utilizing 25% of CPU for Dual core i3 370M processor

2012-08-17 Thread Shantanu MULLICK
Hello Everyone

I have a dual-core Intel i3-370M processor and Windows 64 operating system.
I have a 3GB RAM and my processor can support a maximum of 8GB RAM.

I have virtual memory enabled on my computer.

I am running a program in the 64 bit R which implements a MCMC on a large
dataset and involves around 8 iterations.

The processor estimates that it will take around 1000 minutes to complete
running the program.

However when I check the CPU usage it shows only 25% usage.

I have read on threads that R runs on only 1 core - and thus was expecting
a CPU usage of around 50%.

It would be great if anyone can point out a way by which I can make the
processor use 50% of the CPU which I believe would allow the program to be
run quicker.

Best
Shantanu

[[alternative HTML version deleted]]

__
R-help@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.


[R] Parameterization of Inverse Wishart distribution available in MCMCpack and bayesm libraries

2012-03-01 Thread Shantanu MULLICK
Hello Everyone

Both the MCMCpack and the bayesm libraries allow us to make draws from the
Inverse Wishart distribution.

But I wanted to find out how exactly is the Inverse Wishart distribution
parameterized in these libraries.

The reason I ask is the following:

Now its generally standard to express Inverse Wishart  as IW(0.5 * DOF,0.5*
Scale).(DOF- Degree of freedom, Scale - Scale parameter).

If we follow standard usage when we refer to the Degree of Freedom of the
above IW distribution it is =  DOF (and not 0.5* DOF).

Similarly the Scale parameters of the above IW it is= Scale (and not
0.5*Scale).

For the MCMCpack the IW draws are made by riwish(v,S).

*Question:* Does this compute IW( v,S) or IW(0.5*v,0.5*S) ?

This is the reason I want to find out the way these libraries parameterize
the Inverted Wishart distribution.

Best
Shantanu

[[alternative HTML version deleted]]

__
R-help@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.


[R] System is computationally singular error when using cholesky decompostion in MCMC

2012-02-20 Thread Shantanu MULLICK
Hello Everyone

I have a MCMC loop to calculate a time varying hierarchical Bayesian
structure.

This requires me to use around 5-6 matrix inversions in the loop.

I use cholesky and chol2inv for the matrix decomposition.

Because of the data I am working with I am required to invert a 167 by 167
matrix twice in one iteration.

I need to run the iteration for 1 times, but I get the error  System
is computationally singular  after 5 iterations of the MCMC.

Are there any usual suspects, frequently occurring errors that I should
pay attention to ? My code is 200 lines so I don't write it here.

Naive question 1: Is it possible that the error is occurring for any other
reason other than matrix inversion ?

Naive question 2: Any suggestions on which debugging function (out of the
many functions that R provides)  would be the best to detect the error
given my programming setup. I haven't used any of the debugging functions
in R yet - in the past have been able to debug by going through the code.

Thanks a lot !

Best
Shantanu

[[alternative HTML version deleted]]

__
R-help@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.


[R] Simulating from a Normal Inverted Wishart distribution

2012-02-05 Thread Shantanu MULLICK
Hello everyone

I was wondering how would one simulate from a Normal Wishart Distribution
in R.

A normal inverted Wishart distribution is denoted by NIW (M,C,d,S), where

X/(Sigma) ~ N(  M,C,(Sigma) ) - a matrix normal distribution, (Sigma) -
column dispersion matrix

(Sigma) ~ IW (d,S) - inverted Wishart distribution

Thanks a lot !

Best
Shantanu

[[alternative HTML version deleted]]

__
R-help@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.


Re: [R] Simulating from matrix variate normal distribution

2012-02-05 Thread Shantanu MULLICK
This was indeed very helpful - thanks a lot ! :-))


On 4 February 2012 14:56, Ranjan Maitra mai...@iastate.edu wrote:

 I think this should be easy to write a function doing this.

 Assume that Y is matrix normal with mean matrix mu and row and column
 dispersion matrices Sigma and Gamma, respectively.

 Isn't Y = AZB + mu, where Z is a matrix of independent N(0, 1)'s, A is
 the square root matrix of Sigma (the dispersion matrix of the rows) and
 B is the square root matrix of Gamma. It should be easy to write this
 function in R.

 Again, the density should be easy to write. But it is not clear what
 you mean by a quantile function (of a matrix variate distribution). A
 cdf is going to be a lot harder, but still doable.

 Hope this helps!

 best wishes,
 Ranjan





 On Sat, 4 Feb 2012 00:57:45 +0100 Shantanu MULLICK
 b00295...@essec.edu wrote:

  Hello everyone
 
  Is there a function/command to simulate from matrix variate normal
  distribution in R.
 
  A follow up question would be is there a function/command to obtain the
  density, distribution and quantile function of matrix variate normal
  distribution in R.
 
  Wikipedia has a good description of matrix variate normal distribution
  which is also alternatively called matrix normal distribution.
 
  Thanks a lot !
 
  Best
  Shantanu
 
[[alternative HTML version deleted]]
 
  __
  R-help@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.
 

 __
 R-help@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.


[[alternative HTML version deleted]]

__
R-help@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.


[R] Simulating from matrix variate normal distribution

2012-02-03 Thread Shantanu MULLICK
Hello everyone

Is there a function/command to simulate from matrix variate normal
distribution in R.

A follow up question would be is there a function/command to obtain the
density, distribution and quantile function of matrix variate normal
distribution in R.

Wikipedia has a good description of matrix variate normal distribution
which is also alternatively called matrix normal distribution.

Thanks a lot !

Best
Shantanu

[[alternative HTML version deleted]]

__
R-help@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.