[Rd] [R-SIG-Mac] rnorm.halton

Sat Oct 10 18:04:48 CEST 2009

Hi all,

I need to transform classic 32bit Fortran code to 64bit Fortran code,  
see the discussion [R-SIG-Mac] rnorm.halton. But I'm clearly a  
beginner in Fortran...

Does someone already do this for his package?

 From here, http://techpubs.sgi.com/library/tpl/cgi-bin/getdoc.cgi?coll=linux&db=bks&fname=/SGI_Developer/Porting_Guide/ch03.html 
  , I identify the following changes (Fortran types) to the Fortran  
code:

   - INTEGER to INTEGER*8
   - REAL*8 to REAL*16

The code I would like to change is available on R forge here : http://r-forge.r-project.org/plugins/scmsvn/viewcvs.php/pkg/randtoolbox/src/LowDiscrepancy.f?rev=4229&root=rmetrics&view=markup

Another question I have is how do I tell the R code to use the 64bit  
version of my code on 64bit machine?

In the current implementation, the file quasiRNG.R calls directly the  
Fortran code with .Fortran.
How could I use the 64bit version directly in the R code?

I suspect I need to have a quasiRNG.c file where I will use  
preprocessor statement that will select the good version of the  
function to call. Is that correct?

Thanks in advance

Christophe

Le 15 sept. 2009 à 18:25, Anirban Mukherjee a écrit :

> Sorry, what I should have said was Halton numbers are quasi-random,
> and not pseudo-random. Quasi-random is the technically appropriate
> terminology.
>
> Halton sequences are low discrepancy: the subsequence looks/smells
> random. Hence, they are often used in quasi monte carlo simulations.
> To be precise, there is only 1 Halton sequence for a particular prime.
> Repeated calls to Halton should return the same numbers. The first
> column is the Halton sequence for 2. the second for 3, the third for 5
> and so on using the first M primes (for M columns). (You can also
> scramble the sequence to avoid this.)
>
> I am using them to integrate over a multivariate normal space. If you
> take 1000 random draws, then sum f() over the draws is the expectation
> of f(). If f() is very non-linear (and/or multi-variate) then even
> with large N, its often hard to get a good integral. With quasi-random
> draws, the integration is better for the same N. (One uses the inverse
> distribution function.) For an example, you can look at Train's paper
> (page 4 and 5 have a good explanation) at:
>
> http://elsa.berkeley.edu/wp/train0899.pdf
>
> In the context of simulated maximum likelihood estimation, such
> integrals are very common. Of course true randomness has its own
> place/importance: its just that quasi-random numbers can be very
> useful in certain contexts.
>
> Regards,
> Anirban
>
> PS: qnorm(halton()) gets around the problem of the random deviates  
> not working.
>
> On Tue, Sep 15, 2009 at 11:37 AM, David Winsemius
> <dwinsemius at comcast.net> wrote:
>>
>> On Sep 15, 2009, at 11:10 AM, Anirban Mukherjee wrote:
>>
>>> Thanks everyone for your replies. Particularly David.
>>>
>>> The numbers are pseudo-random. Repeated calls should/would give the
>>> same output.
>>
>> As I said, this package is not one with which I have experience. It
>> has _not_ however the case that repeated calls to (typical?) random
>> number functions give the same output when called repeatedly:
>>
>>  > rnorm(10)
>>   [1] -0.8740195  2.1827411 -0.1473012 -1.4406262  0.1820631
>> -1.3151244 -0.4813703  0.8177692
>>   [9]  0.2076117  1.8697418
>>  > rnorm(10)
>>   [1] -0.7725731  0.8696742 -0.4907099  0.1561859  0.5913528
>> -0.8441891  0.2285653 -0.1231755
>>   [9]  0.5190459 -0.7803617
>>  > rnorm(10)
>>   [1] -0.9585881 -0.0458582  1.1967342  0.6421980 -0.5290280
>> -1.0735112  0.6346301  0.2685760
>>   [9]  1.5767800  1.0864515
>>  > rnorm(10)
>>   [1] -0.60400852 -0.06611533  1.00787048  1.48289305  0.54658888
>> -0.67630052  0.52664127 -0.36449997
>>   [9]  0.88039397  0.56929333
>>
>> I cannot imagine a situation where one would _want_ the output to be
>> the same on repeated calls unless one reset a seed. Unless perhaps I
>> am not understanding the meaning of "random" in the financial domain?
>>
>> --
>> David
>>
>>>   Currently, Halton works fine when used to just get the
>>> Halton sequence, but the random deviates call is not working in 64  
>>> bit
>>> R. For now, I will generate the numbers in 32 bit R, save them and
>>> then load them back in when using 64 bit R. The package maintainers
>>> can look at it if/when they get a chance and/or access to 64 bit R.
>>>
>>> Thanks!
>>>
>>> Best,
>>> Anirban
>>>
>>> On Tue, Sep 15, 2009 at 9:01 AM, David Winsemius <dwinsemius at comcast.net
>>>> wrote:
>>>> I get very different output from the two versions of Mac OSX R as
>>>> well. The 32 bit version puts out a histogram that has an expected,
>>>> almost symmetric unimodal distribution. The 64 bit version  
>>>> created a
>>>> bimodal distribution with one large mode near 0 and another smaller
>>>> mode near 10E+37. Postcript output attached.
>>>>
>>>
>>>
>>>
>>> --
>>> Anirban Mukherjee | Assistant Professor, Marketing | LKCSB, SMU
>>> 5062 School of Business, 50 Stamford Road, Singapore 178899 |
>>> +65-6828-1932
>>>
>>> _______________________________________________
>>> R-SIG-Mac mailing list
>>> R-SIG-Mac at stat.math.ethz.ch
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mac
>>
>> David Winsemius, MD
>> Heritage Laboratories
>> West Hartford, CT
>>
>>
>
>
>
> -- 
> Anirban Mukherjee | Assistant Professor, Marketing | LKCSB, SMU
> 5062 School of Business, 50 Stamford Road, Singapore 178899 |  
> +65-6828-1932
>
> _______________________________________________
> R-SIG-Mac mailing list
> R-SIG-Mac at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-mac

--
Christophe Dutang
Ph.D. student at ISFA, Lyon, France
website: http://dutangc.free.fr