How to solve a expression with 2 variables that are within pchisq() in R?

rstudio
statistics
functions

#1

Hi everyone,

I am the newcomer and trying to solve a math expression that have two variables, ca and y , and actually I only want to know the value of variable ca because the variable y will be integrated like the code below:

#---time start---#
ptm <- proc.time()

library(pracma)

U = 14; L = 7; d = (U-L)/2;
T = ( 3*U+L )/(3+1); #  by 3Du = Dl 
Du = U-T; Dl = T - L; d_star = min(Du,Dl);
du = d/Du; dl = d/Dl;
alpha = 0.05;
aLe = 0.05; 
n = 100;

ksi_given = -3;

if (ksi_given <= 0){
  sigma = fzero(function(x) aLe - (x^2)*(  (dl^2*ksi_given^2+1)/d_star^2  ), c(0,0.5) )$x ;
  mu = ksi_given*sigma + T;
  B = n*(  (dl*ksi_given)^2 + 1  )/aLe;
  C = (dl*ksi_given)^2/(d_star/sigma)^2 + 1/(d_star/sigma)^2;
  b = d_star/sigma;
}else{
  sigma = fzero(function(x) aLe - (x^2)*(  (du^2*ksi_given^2+1)/d_star^2  ), c(0,0.5) )$x ;
  mu = ksi_given*sigma + T;
  B = (   n*(  (du*ksi_given)^2 + 1  )   )/aLe;
  C = (du*ksi_given)^2/(d_star/sigma)^2 + 1/(d_star/sigma)^2;
  b = d_star/sigma;
}

fun = function(ca,y) { pchisq( B*ca-y, n-1)* ( 1/ (2*sqrt(y)) )* (   (dl^-1)*normpdf( (dl^-1)*sqrt(y) + sqrt(n*ksi_given) ) + (du^-1)*normpdf( (du^-1)*sqrt(y) - sqrt(n*ksi_given) )    ) };

fun2 = function(ca)  alpha - integrate(fun(ca,y),lower = 0.0001,upper = B*ca);

c = linspace(1e-4,0.1,100);
for (i in 1:numel(c)){
  fc(i) = fun2(c(i));
}

plot(c,fc)

ca = fzero(function(ca) fun2, aLe);

About the expression that I will use fzero() to solve, it is created by these steps:

  1. There is an originally expression fun which has two variables, ca and y.
  2. I will use integrate() to integrate fun about the variable y, and then get the fun2 .

Finally, I will try to use fzero() to solve the value of ca in fun2 like my last line of code.

Unfortunately, R return the error message like this:

Error in pchisq(B * ca - y, n - 1) : argument "y" is missing, with no default

I have also tried this modifying code

ca = fzero(function(ca) alpha - integrate({ pchisq( B*ca-y, n-1)* ( 1/ (2*sqrt(y)) )* (   (dl^-1)*normpdf( (dl^-1)*sqrt(y) + sqrt(n)*(ksi_given) ) + (du^-1)*normpdf( (du^-1)*sqrt(y) - sqrt(n)*(ksi_given) )    ) },lower = 0.0001,upper = B*ca), aLe);

but it doesn't work.

Does anyone can give me some advice about my problem to correctly get the value of ca?

Does it mean pchisq() could not using variables as its parameters? or I have misunderstand the usage of function() or integrate() or fzero()?

Hope you can give me some advice! Thanks!


#2

This does not look like something I'd use R for. From time to time, I too need to solve equations, usually I turn to pen-on-paper, but if it's too complicated and time consuming, I use wolframalpha

If you click the link at write solve(a + b*x + c*x^2, x), then it'll do just that for you.

Hope it helps! :slightly_smiling_face:


#3

Hi Leon,

You are right. It seems like the function to solve a equation is not the main function in R and actually, this equation it is so hard to solve by pen-on-paper method, did you have any idea?

Thanks!


#4

As I suggested, take a look at Wolfram Alpha?


#5

I have tried it, but the result is not what I expected.


#6

What did you enter? You can integrate like so integrate(1/x, 2, 4) will integrate 1/x from lower limit 2, to upper limit 4, which is equal to `log(2)


#7

Hi Leon,

The original expression that I want to solve is so complicated like this:

It is not so normal like the example that you give...


#8

Depending on Phi and F_K, it's not too bad - I'd do pen on paper to break the expression down into a series of additions and then use WolframAlpha to solve each part of the expression :slightly_smiling_face:

But either way algebra like this is definitely not a job for R

Hope it helps!


#9

Thanks! I will try other methods to solve the algebra problem!