2

I'm trying to compute the natural log of the exponential plus the median of an array but precision of the floating point numbers has to be at 2000 or the answer will always be 0.

Here is what I have so far:

import bigfloat    
x = np.array([-15349.79, -15266.66, -15242.86])
answer = np.median(x) - bigfloat.log(np.mean(bigfloat.exp(x, precision(2000)) + np.median(x)))

This code returns the following error because BigFloat.exp doesn't work on list types.

__main__:1: RuntimeWarning: invalid value encountered in log
    array([ nan,  nan,  nan])
    >>> np.median - bigfloat.log(np.mean(bigfloat.exp(x, precision(2000)) + np.median(x)))
    Traceback (most recent call last):
      File "<stdin>", line 1, in <module>
      File "/usr/local/lib/python2.7/dist-packages/bigfloat/core.py", line 1446, in exp
        (BigFloat._implicit_convert(x),),
      File "/usr/local/lib/python2.7/dist-packages/bigfloat/core.py", line 800, in _implicit_convert
        "to BigFloat" % (arg, type(arg)))
    TypeError: Unable to convert argument [-15349.79 -15266.66 -15242.86] of type <type 'numpy.ndarray'> to BigFloat

I then tried list comprehension to build exponentials with precision=2000 but the output loses precision.

>>>[bigfloat.exp(num, precision(2000)) + BigFloat(np.median(x), context=precision(2000)) for num in x]

[BigFloat.exact('-15266.660000000000', precision=53), BigFloat.exact('-15266.660000000000', precision=53), BigFloat.exact('-15266.660000000000', precision=53)]

I can get exponential values from x with precision using list comprehension like this:

>>> y = np.array([bigfloat.exp(num, precision(2000)) for num in x])
array([ BigFloat.exact('4.687104391719151845495683439738733373338442839115295629901427982078731632624437676931490944129698360820181432182432247030174990996132679553420878801399467284014932435857939862423137788809352433385343338286781879066456873472082636971652587852379587642184738259777920214145750081457830886095059600271300801524086458081526778407096875577469007859748418787409261916565428551066536598870494431653463691369520685259389041691328858076222360659135651318185915247243253437425313718584237418827954902711646850362440592578566737383074556346759332830833101196990845567845361685702120838400980658468192376607029699380e-6667', precision=2000),
       BigFloat.exact('5.940252515613784074556309936460240719247289162721506036367173706446251980105141516508485353816420203092949808160441872102510858074752258842471232457791645630942356791030979819957875784926363292788234347364203226596687259399429746917920741954706671428927952815041593801046941161417184730521922391569252658359604561007344640221209247257690437938341582528542446104212627600044421405229891732857592357725820529732658234033631318425476205945557122313924830225909708184718714804450362854280687393958953216361511561704836206669636791590655443171684811683402330053964681285971765503804198242627954802184910024902e-6631', precision=2000),
       BigFloat.exact('1.288289823457680769007768910906908351089465066737359292851136781233384628626343352992636825492031571595103435334370758196286825224816434202324210124540257467624499933926757723584914581286627071375803686797647455160335628799775858142489656885909172399293492295645703203222218182084249732700966821340431452575815746282651823925677898549539269454644804048695031225101971491985645951419388454738446335860070391719831039721150295437452237572831818257568221336263814212095143032293694959570063677473370798577031762607527221992020661837810811705007892502559603984057933079724597470162879821292957122400421187875e-6620', precision=2000)], dtype=object)

So, how do I add np.median(x) to each item in the exponential array and how do I get the log of the elements in the final array? Is there an easier way to compute the equation given by variable answer in the first code snippet above?

I'm basically trying to convert this R code to Python:

llMed <- stats::median(x)
metric <- as.double(
  llMed - log( Rmpfr::mean( exp( -Rmpfr::mpfr(x, prec=2000L) + llMed )))
)
Kenneth Orton
  • 399
  • 1
  • 11

1 Answers1

0

I think I was able to transpose the R code above, using gmpy2 instead of BigFloat but probably could have used either:

import gmpy2 as gp
import numpy as np    
with gp.local_context(gp.context(), precision=2000) as ctx:
    x = np.array([gp.mpfr(-15349.79), gp.mpfr(-15266.66), gp.mpfr(-15242.86)])
    answer = np.median(x) - gp.log(np.mean([gp.exp(gp.add(np.median(x), -gp.mpfr(item))) for item in x]))

Output returns:

mpfr('-15348.69138771133276342351845677418587273364155071266454924792376077085597654860712492594599688058199377288639137666410888137048513058096745776113277020531535761588878042016412651412489608285582998650481415959482636259292554205272401992800167437149224493900214813760731711328834885525224938584036572786202764947113186177391441168267495103157033223800214492991771147327262071249020907408433551462034872024117419365924381891832161483692689444158313503155805562703358104122358438183824459422779020196496507161021965435781332319270106503099589290788685588200038739438059696970511568363022088057740627040541967',2000)
Kenneth Orton
  • 399
  • 1
  • 11