CODEWITHSUNDEEP
r
user2884679
user2884679

Reputation: 129

Weird error in R

For some reason, the if-statement in my code generates some weird error in R whenever I run my function: 

Error in if (P > 0) { : missing value where TRUE/FALSE needed

My function is supposed to return the smallest eigenvalue of the symmetric, positive semi-definite matrix A, if B is the identity matrix. I tested my function using

set.seed(12345)
a <- crossprod(matrix (rnorm(40), 10 ,4)) / 10
b <- diag(1, 4, 4)

and it does yield the correct answer (I checked using the eigen function.)

Anyway, I need to set theta to the value that minimizes a certain function. For the example I used, that value is (-Q - sqrt(under.sqrt)) / (2 * P), but in general it depends on the sign of P. Why does my if-statement give me such an error? I have been stuck here for a while. Any help would be appreciated. 

myFunction <- function(A, B, eps = 1e-6, itmax = 100, verbose = FALSE) {

  n <- nrow(A)
  m <- ncol(A)
  x <- rep (1, n)
  u <- A %*% x
  v <- B %*% x
  r <- t(u) %*% x
  s <- t(v) %*% x
  y <- r / s

  itel <- 1

  repeat {
    tmax <- 0

    for(j in 1:m) {
      P <- (u[j] * B[j, j]) - (v[j] * A[j, j])
      Q <- (r * B[j, j]) - (s * A[j, j])
      R <- (r * v[j]) - (s * u[j])

      under.sqrt <- Q^2 - 4 * P * R

      # Error right here
      if (P > 0) {
        theta <- (-Q + sqrt(under.sqrt)) / (2 * P)
      }

      else {
        theta <- (-Q - sqrt(under.sqrt)) / (2 * P)
      }


      x[j] <- x[j] + theta
      r <- r + 2 * theta * u[j] + A[j, j] * theta * theta
      s <- s + 2 * theta * v[j] + B[j, j] * theta * theta
      u <- u + theta * A[ ,j]
      v <- v + theta * B[, j]

      y <- r / s

      tmax <- max(tmax, abs(theta))
    }

      if ((tmax < eps) || (itel == itmax)) {
        break()
      }

      itel <- itel + 1  

  }

  return(y)
}

Upvotes: 2

Views: 2737

Answers (4)

James
James

Reputation: 66834

You algorithm is not stable. Using options(error=recover) you can browse the function where it is going wrong:

Browse[1]> P
[1] NaN
Browse[1]> theta
     [,1]
[1,]  Inf

Here you can see that P is NaN and this is caused by theta blowing up to infinite values. This feeds through to x, u, v and ultimately P and makes the comparison fail.

Upvotes: 4

Raffael
Raffael

Reputation: 20045

Here's my bet: P is NA

> if(NA > 0){}
Error in if (NA > 0) { : missing value where TRUE/FALSE needed

P renders NaN:

> if(NaN > 0){}
Error in if (NaN > 0) { : missing value where TRUE/FALSE needed

regarding the guess that P is a vector:

> if(1:2 > 0){}
NULL
Warning message:
In if (1:2 > 0) { :
  the condition has length > 1 and only the first element will be used

Upvotes: 2

zx8754
zx8754

Reputation: 56149

Try ifelse:

theta <- ifelse(P > 0,
                (-Q + sqrt(under.sqrt)) / (2 * P),
                (-Q - sqrt(under.sqrt)) / (2 * P))

ifelse can handle NA comparison:

ifelse(NA>0,print("x"),print("y"))
[1] NA

if(NA>0) print("x") else print("y")
Error in if (NA > 0) print("x") else print("y") : 
  missing value where TRUE/FALSE needed

Upvotes: 0

Hans Roelofsen
Hans Roelofsen

Reputation: 741

Are you sure P is a single value? Ie, not a vector of length > 1?

Upvotes: 0

Related Questions