How to Calculate Any Negative Number to the Power of Some Fraction in R

How to calculate any negative number to the power of some fraction in R?

As documented in help("^"):

Users are sometimes surprised by the value returned, for example
why ‘(-8)^(1/3)’ is ‘NaN’. For double inputs, R makes use of IEC
60559 arithmetic on all platforms, together with the C system
function ‘pow’ for the ‘^’ operator. The relevant standards
define the result in many corner cases. In particular, the result
in the example above is mandated by the C99 standard. On many
Unix-alike systems the command ‘man pow’ gives details of the
values in a large number of corner cases.

So you need to do the operations separately:

R> ((-8/27)^2)^(1/3)
[1] 0.4444444

Raising vector with negative numbers to a fractional exponent in R

exponent <- function(a, pow) (abs(a)^pow)*sign(a)

how to generate a sequences of positive and negative fractional numbers in R

The answer below adds the new requirement at the end (about not having 0).

You need seq and to have the numbers in the correct order:

sample(seq(0.2,from=-0.9,by=.1), 10, replace=T)

[1] -0.1 -0.6  0.2  0.0  0.0  0.1  0.2 -0.8 -0.7 -0.2

and I recommend sample_n:

library(dplyr)
sample_n(as_tibble(seq(0.2,from=-0.9,by=.1)), 10, replace=T)

    value
    <dbl>
 1  0.100
 2 -0.200
 3  0.   
 4 -0.100
 5  0.   
 6 -0.100
 7 -0.600
 8 -0.800
 9 -0.900
10 -0.900

Update:

You mentioned you don't want 0. There are many ways to accomplish this. The best way probably depends on the set of statistical implications that you prefer if that's important to you. Here are some examples:

Example 1

x <- sample(seq(0.2,from=-0.9,by=.1), 10, replace=T)

ifelse(x==0, sample(seq(-.1,from=-0.9,by=.1), 1, replace=T), x)

Example 2

c(sample(seq(-.1,from=-0.9,by=.1), 5, replace=T),
  sample(seq(.1, 0.2,by=.1), 5, replace=T))

You could also use a while. Basically for each element you'd take random samples until there's one that's not == 0. That's probably how I would do it to avoid polluting the distributional assumptions.

R expression results in NaN for no obvious reason

It is because of the implementation of pow under C99 standard.

Let alone OP's example: (-50.61828)^(-67.421587), the mathematically justified (-8)^(1/3) = -2 does not work in R:

(-8)^(1/3)
# [1] NaN

Quoted from ?"^":

 Users are sometimes surprised by the value returned, for example
 why ‘(-8)^(1/3)’ is ‘NaN’.  For double inputs, R makes use of IEC
 60559 arithmetic on all platforms, together with the C system
 function ‘pow’ for the ‘^’ operator.  The relevant standards
 define the result in many corner cases.  In particular, the result
 in the example above is mandated by the C99 standard.  On many
 Unix-alike systems the command ‘man pow’ gives details of the
 values in a large number of corner cases.

I am on Ubuntu LINUX, so can help get relevant part of man power printed here:

   If x is a finite value less than 0, and y is  a  finite  noninteger,  a
   domain error occurs, and a NaN is returned.

Rule to calculate power of a number when the exponent is Negative in Prolog?

First, one should consider how to define 0⁰. Formally speaking it is indeterminate. It could be zero or it could be 1. As Wolfram's Mathworld says in its article on powers and in its article on zero:

0⁰ (zero to the zeroth power) itself is undefined. The lack of a well-defined meaning for this quantity follows from the mutually contradictory facts that a⁰ is always 1, so 0⁰ should equal 1, but 0^a is always 0 (for a > 0), so 0^a should equal 0. The choice of definition for 0⁰ is usually defined to be indeterminate, although defining 0⁰ = 1 allows some formulas to be expressed simply (Knuth 1992; Knuth 1997, p. 57).

So you should first choose how to define the special case of 0⁰: Is it 0? Is it 1? Is it undefined?

I choose to look at it as being undefined.

That being said, you can look at a positive exponent as indicated repeated multiplication (e.g. 10³ is 10*10*10, or 1,000), and you can look at a negative exponent as indicating repeated division (e.g, 10^-3 is (((1/10)/10)/10), or 0.001). My inclination, partly because I like the symmetry of this approach and partly to avoid the cuts (since a cut is often a signal that you've not defined the solution properly), would be something like this:

% -----------------------------
% The external/public predicate
% -----------------------------
pow( 0 , 0 , _ ) :- ! , fail .
pow( X , N , R ) :-
  pow( X , N , 1 , R )
  .

% -----------------------------------
% the tail-recursive worker predicate
% -----------------------------------
pow( _ , 0 , R , R  ).
pow( X , N , T , R  ) :-
  N > 0 ,
  T1 is T * X ,
  N1 is N-1   ,
  pow( X , N1 , T1 , R )
  .
pow( _ , 0 , R , R  ) :-
  N < 0 ,
  T1 is T / X ,
  N1 is N+1   ,
  pow( X , N1 , T1 , R )
  .

The other approach, as others have noted, is to define a positive exponent as indicating repeated multiplication, and a negative exponent as indicating the reciprocal of the positive exponent, so 10³ is 10*10*10 or 1,000, and 10^-3 is 1/(10³), or 1/1,000 or 0.001. To use this definition, I'd again avoid the cuts and do something like this:

% -----------------------------
% the external/public predicate
% -----------------------------
pow( 0 , 0 , _ ) :-  % 0^0 is indeterminate. Is it 1? Is it 0? Could be either.
  ! ,
  fail
  .
pow( X , N , R ) :-
  N > 0 ,
  pow( X , N , 1 , R )
  .
pow( X , N , R ) :-
  N < 0 ,
  N1 = - N ,
  pow( X , N1 , 1 , R1 ) ,
  R is 1 / R1
  .

% -----------------------------------
% The tail-recursive worker predicate
% -----------------------------------
pow( _ , 0 , R , R  ).
pow( X , N , T , R  ) :-
  N > 0 ,
  T1 is T * X ,
  N1 is N-1   ,
  pow( X , N1 , T1 , R )
  .

In Python 3.6, why does a negative number to the power of a fraction return nan when in a numpy array?

Exponentiation in python has higher precedence than the negative operator. Thus -1000**(1/3) is equivalent to -(1000**(1/3)).

When you doing this operation inside the loop you get (-1000)**(1/3). This equal to 10 * (-1**(1/3)) which a complex number. Now the array you have, uses a default data type since you did not define any that is determined according to the documentation as follows:

dtype : data-type, optional

The desired data-type for the array. If not given, then the type will be determined as the minimum type required to hold the objects in the sequence. This argument can only be used to ‘upcast’ the array. For downcasting, use the .astype(t) method.

So it is probably np.int16.

Putting all the information together, we can conclude that your array is not equipped with the appropriate dtype attribute to be able to hold the result of (-1000)**(1/3) even though the result exists.

This does not happen outside arrays since there, no dtype is assumed.

---

Fix \ Workaround:

>>> a = np.array([-1000, 1], dtype=np.complex)
>>> for i in a:
...     print(i**(1/3.))
...
(5+8.66025403784j)
(1+0j)

---

Related Topics