Store the output of partial derivative

timtimalo · January 5, 2018, 1:16am

I calculated partial derivative of functions F2 and F6 with respect to unknown; a and b and I want to store the results of partial derivative in a matrix called “mat” to use it later in the R-program.
I tried the following code to calculate the partial derivatives and store the results in a matrix “mat” that has 2 columns and 2 rows. I got an error says:

Error: object 'a' not found
Error: object 'b' not found

So; what is the suitable code to store the results of the derivative of functions with respect to a and b in matrix mat and then do some matrix work using matrix mat? I appreciate your help.
The codes

x=c("a","b")#vector of unknown parameters; a and b
r=matrix(0,nrow=2,ncol=2)#create matrix 2 by 2
library("Deriv")
F2=function(a,b) (choose(10,0)*beta((a*(1-b))/b+0,(1-a)*(1-b)/b+10)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                    choose(10,1)*beta((a*(1-b))/b+1,(1-a)*(1-b)/b+9)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                    choose(10,2)*beta((a*(1-b))/b+2,(1-a)*(1-b)/b+8)/beta((a*(1-b))/b,(1-a)*(1-b)/b))
Deriv(F2,x)
mat[1,1]=c(a)
mat[1,2]=c(b)
F6=function (a,b) (choose(10,0)*beta((a*(1-b))/b+0,(1-a)*(1-b)/b+10)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,1)*beta((a*(1-b))/b+1,(1-a)*(1-b)/b+9)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,2)*beta((a*(1-b))/b+2,(1-a)*(1-b)/b+8)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,3)*beta((a*(1-b))/b+3,(1-a)*(1-b)/b+7)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,4)*beta((a*(1-b))/b+4,(1-a)*(1-b)/b+6)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,5)*beta((a*(1-b))/b+5,(1-a)*(1-b)/b+5)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,6)*beta((a*(1-b))/b+6,(1-a)*(1-b)/b+4)/beta((a*(1-b))/b,(1-a)*(1-b)/b))
Deriv(F6,x)
mat[2,1]=c(a)
mat[2,2]=c(b)
mat

The results

function (a, b) 
{
    .e1 <- 1 - b
    .e2 <- 1 - a
    .e4 <- .e2 * .e1/b
    .e6 <- a * .e1/b
    .e7 <- .e1/b
    .e9 <- digamma(.e4)
    .e10 <- digamma(.e7 + 10)
    .e11 <- digamma(.e7)
    .e12 <- .e4 + 10
    .e13 <- .e4 + 8
    .e14 <- .e4 + 9
    .e15 <- 1 + .e6
    .e16 <- 2 + .e6
    .e17 <- digamma(.e6)
    .e19 <- .e2 * (.e9 - .e11) + a * (.e17 - .e11)
    .e20 <- b * beta(.e6, .e4)
    .e21 <- beta(.e15, .e14)
    .e22 <- beta(.e16, .e13)
    .e23 <- beta(.e6, .e12)
    .e24 <- digamma(.e12)
    .e25 <- digamma(.e13)
    .e26 <- digamma(.e14)
    .e27 <- digamma(.e15)
    .e28 <- digamma(.e16)
    .e29 <- .e17 - .e9
    c(a = ((10 * (.e27 - .e26) - 10 * .e29) * .e21 + (45 * (.e28 - 
        .e25) - 45 * .e29) * .e22 + .e23 * (.e9 - .e24)) * .e1/.e20, 
        b = ((.e2 * (.e9 + .e10 - (.e24 + .e11)) + a * (.e10 - 
            .e11)) * .e23 + (10 * .e19 - 10 * (.e2 * (.e26 - 
            .e10) + a * (.e27 - .e10))) * .e21 + (45 * .e19 - 
            45 * (.e2 * (.e25 - .e10) + a * (.e28 - .e10))) * 
            .e22) * (.e7 + 1)/.e20)
}

mat[1,1]=c(a)
# Error: object 'a' not found
mat[1,2]=c(b)
# Error: object 'b' not found
F6=function (a,b) (choose(10,0)*beta((a*(1-b))/b+0,(1-a)*(1-b)/b+10)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,1)*beta((a*(1-b))/b+1,(1-a)*(1-b)/b+9)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,2)*beta((a*(1-b))/b+2,(1-a)*(1-b)/b+8)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,3)*beta((a*(1-b))/b+3,(1-a)*(1-b)/b+7)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,4)*beta((a*(1-b))/b+4,(1-a)*(1-b)/b+6)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,5)*beta((a*(1-b))/b+5,(1-a)*(1-b)/b+5)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,6)*beta((a*(1-b))/b+6,(1-a)*(1-b)/b+4)/beta((a*(1-b))/b,(1-a)*(1-b)/b))
Deriv(F6,x)
function (a, b) 
{
    .e1 <- 1 - b
    .e2 <- 1 - a
    .e4 <- .e2 * .e1/b
    .e6 <- a * .e1/b
    .e7 <- .e1/b
    .e9 <- digamma(.e4)
    .e10 <- digamma(.e7 + 10)
    .e11 <- digamma(.e7)
    .e12 <- digamma(.e6)
    .e14 <- .e2 * (.e9 - .e11) + a * (.e12 - .e11)
    .e15 <- .e12 - .e9
    .e16 <- .e4 + 10
    .e17 <- .e4 + 4
    .e18 <- .e4 + 5
    .e19 <- .e4 + 6
    .e20 <- .e4 + 7
    .e21 <- .e4 + 8
    .e22 <- .e4 + 9
    .e23 <- 1 + .e6
    .e24 <- 2 + .e6
    .e25 <- 3 + .e6
    .e26 <- 4 + .e6
    .e27 <- 5 + .e6
    .e28 <- 6 + .e6
    .e29 <- 210 * .e14
    .e30 <- 210 * .e15
    .e31 <- b * beta(.e6, .e4)
    .e32 <- beta(.e23, .e22)
    .e33 <- beta(.e24, .e21)
    .e34 <- beta(.e25, .e20)
    .e35 <- beta(.e26, .e19)
    .e36 <- beta(.e27, .e18)
    .e37 <- beta(.e28, .e17)
    .e38 <- beta(.e6, .e16)
    .e39 <- digamma(.e16)
    .e40 <- digamma(.e17)
    .e41 <- digamma(.e18)
    .e42 <- digamma(.e19)
    .e43 <- digamma(.e20)
    .e44 <- digamma(.e21)
    .e45 <- digamma(.e22)
    .e46 <- digamma(.e23)
    .e47 <- digamma(.e24)
    .e48 <- digamma(.e25)
    .e49 <- digamma(.e26)
    .e50 <- digamma(.e27)
    .e51 <- digamma(.e28)
    c(a = ((10 * (.e46 - .e45) - 10 * .e15) * .e32 + (120 * (.e48 - 
        .e43) - 120 * .e15) * .e34 + (210 * (.e49 - .e42) - .e30) * 
        .e35 + (210 * (.e51 - .e40) - .e30) * .e37 + (252 * (.e50 - 
        .e41) - 252 * .e15) * .e36 + (45 * (.e47 - .e44) - 45 * 
        .e15) * .e33 + .e38 * (.e9 - .e39)) * .e1/.e31, b = ((.e2 * 
        (.e9 + .e10 - (.e39 + .e11)) + a * (.e10 - .e11)) * .e38 + 
        (10 * .e14 - 10 * (.e2 * (.e45 - .e10) + a * (.e46 - 
            .e10))) * .e32 + (120 * .e14 - 120 * (.e2 * (.e43 - 
        .e10) + a * (.e48 - .e10))) * .e34 + (.e29 - 210 * (.e2 * 
        (.e40 - .e10) + a * (.e51 - .e10))) * .e37 + (.e29 - 
        210 * (.e2 * (.e42 - .e10) + a * (.e49 - .e10))) * .e35 + 
        (252 * .e14 - 252 * (.e2 * (.e41 - .e10) + a * (.e50 - 
            .e10))) * .e36 + (45 * .e14 - 45 * (.e2 * (.e44 - 
        .e10) + a * (.e47 - .e10))) * .e33) * (.e7 + 1)/.e31)
}

mat[2,1]=c(a)
# Error: object 'a' not found
mat[2,2]=c(b)
# Error: object 'b' not found
mat
#      [,1] [,2]
# [1,]    0    0

danr · January 5, 2018, 11:55am

It would help if you posted a reprex of your code instead of copied console output, it would make it a lot easier for us to help you. As is we can't even copy/paste you code without fixing up all the printer quotes.

https://www.tidyverse.org/help/#reprex

danr · January 5, 2018, 1:34pm

I don't see anywhere you have assigned values to a or b, that's why the error is occurring. Are you trying to store some functions with some precomputed values in the matrix?

What are you expecting
choose(10,0)beta(...)
to do?

beta(...) returns a number and choose(10,0) returns 1 so what you have here is the equivalent of:

suppressPackageStartupMessages(library(Deriv))
choose(10,0)beta(1.2, 3)
#> Error: <text>:2:13: unexpected symbol
#> 1: suppressPackageStartupMessages(library(Deriv))
#> 2: choose(10,0)beta
#>                ^

which is a syntax error.

That syntax error is probably being masked by the error that is due to a and b not being assigned a value.

Is this what you wanted?

suppressPackageStartupMessages(library(Deriv))
choose(10,0)*beta(1.2, 3)
#> [1] 0.2367424

Even if your script is patched up like this the problem of a and b being unassigned is still a problem that will cause a runtime error.

timtimalo · January 5, 2018, 5:01pm

“a” and “b” are unknown. I have two cumulative distribution functions and each cumulative distribution function is calculated at specific value of random variable; x and they depend on beta functions; “beta(a,b)” which are functions of unknown parameters “a” and “b”. I calculated the partial derivative of both cumulative distribution functions with respect to “a” and “b” using package “Deriv” and I got the results of the derivatives in the form C(a=…., b=….). I want to store the results of partial derivatives for “a” and “b” in a matrix 2 by 2 to do some matrix calculations. I tried to construct an empty matrix; “mat” to store the results of the partial derivatives in it but I got error. How I can correctly store the results of the partial derivatives of the cumulative distribution functions with respect to “a” and “b” in a matrix? Thank you.

I tried to use "reprex" but I could not continue and this what I did and where I stopped

** I installed reprexas**

install.packages("reprex")    
library(reprex)

I wrote the code and copy it to the clipboard
For this step, I do not know what it is meant by paste it in GitHub issue. How I can do it?
I used Addins and the chose Render repex.
I click Render. For where is repex source, I chose on the clipboard. For the source, I chose GitHub.
I do not know what I should do next to put the code in repex?

Thus, I copied the codes here

The following is the code:

x=c("a","b")
mat=matrix(0,nrow=2,ncol=2)#create 2 by 2 matrix
library("Deriv")
# Cumulative distribution function that is function of unknown "a" and "b" and its partial derivative 
# with respect to "a" and "b"
F2=function(a,b) (choose(10,0)*beta((a*(1-b))/b+0,(1-a)*(1-b)/b+10)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                    choose(10,1)*beta((a*(1-b))/b+1,(1-a)*(1-b)/b+9)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                    choose(10,2)*beta((a*(1-b))/b+2,(1-a)*(1-b)/b+8)/beta((a*(1-b))/b,(1-a)*(1-b)/b))
Deriv(F2,x)
# store the results of partial derivative in matrix; mat in the first row
mat[1,1]=c(a)
mat[1,2]=c(b)
# Another cumulative distribution function that is function of unknown "a" and "b" and its partial 
# derivative with respect to "a" and "b" 
F6=function (a,b) (choose(10,0)*beta((a*(1-b))/b+0,(1-a)*(1-b)/b+10)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,1)*beta((a*(1-b))/b+1,(1-a)*(1-b)/b+9)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,2)*beta((a*(1-b))/b+2,(1-a)*(1-b)/b+8)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,3)*beta((a*(1-b))/b+3,(1-a)*(1-b)/b+7)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,4)*beta((a*(1-b))/b+4,(1-a)*(1-b)/b+6)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,5)*beta((a*(1-b))/b+5,(1-a)*(1-b)/b+5)/beta((a*(1-b))/b,(1-a)*(1-b)/b)+
                     choose(10,6)*beta((a*(1-b))/b+6,(1-a)*(1-b)/b+4)/beta((a*(1-b))/b,(1-a)*(1-b)/b))
Deriv(F6,x)
#store the results of partial derivative in matrix; mat in the second row
mat[2,1]=c(a)
mat[2,2]=c(b)
#Show matrix mat
mat

Thank you

danr · January 5, 2018, 5:16pm

I'm out the door in just a few minutes but me or someone else, I'm sure, will take a closer look at what you are trying to do later.

As far a reprex goes... When I use reprex i select the code I want to reprex in the editor and use "select", not "clipboard" as the source. Then I use github for output. It won't send anything to github it just formats in a way that github and community.rstudio.com can handle.

What reprex is doing, basically, is capturing the code you selected along with the output it sends to the console when executed. Then if converts all that to HTML so that it an be pasted into a web page.

So once reprex completes and you press the OK button that "html" is on you clipboard an you can just paste it into a message here and it will end up nicely formatted and it will be code that is know to be runable.

Later

danr · January 5, 2018, 5:27pm

I acted to quickly... just saw the door opening ... this has some errors in it

Now I'm seconds from going out door but after a quick look... This works but I'm not sure it does what you want (I still need to read closer your previous post)

A added the *s between the choose and beta. There is no variable mat but the one named r looks like what you wanted to use.

More later...

suppressPackageStartupMessages(library(Deriv))
x=c("a","b")#vector of unknown parameters; a and b
r=matrix(0,nrow=2,ncol=2)#create matrix 2 by 2
library("Deriv")
F2=function(a,b) (choose(10,0)*beta((a(1-b))/b+0,(1-a)(1-b)/b+10)/beta((a(1-b))/b,(1-a)(1-b)/b)+
                                        choose(10,1)*beta((a(1-b))/b+1,(1-a)(1-b)/b+9)/beta((a*(1-b))/b,(1-a)(1-b)/b)+
                                        choose(10,2)*beta((a(1-b))/b+2,(1-a)(1-b)/b+8)/beta((a*(1-b))/b,(1-a)(1-b)/b))
Deriv(F2,x)
#> Error in simplifications[[sym.name]]: wrong arguments for subsetting an environment
mat[1,1]=c(a)
#> Error in eval(expr, envir, enclos): object 'a' not found
mat[1,2]=c(b)
#> Error in eval(expr, envir, enclos): object 'b' not found
F6=function (a,b) (choose(10,0)*beta((a(1-b))/b+0,(1-a)(1-b)/b+10)/beta((a*(1-b))/b,(1-a)(1-b)/b)+
                                        choose(10,1)*beta((a(1-b))/b+1,(1-a)(1-b)/b+9)/beta((a*(1-b))/b,(1-a)(1-b)/b)+
                                        choose(10,2)*beta((a(1-b))/b+2,(1-a)(1-b)/b+8)/beta((a*(1-b))/b,(1-a)(1-b)/b)+
                                        choose(10,3)*beta((a(1-b))/b+3,(1-a)(1-b)/b+7)/beta((a*(1-b))/b,(1-a)(1-b)/b)+
                                        choose(10,4)*beta((a(1-b))/b+4,(1-a)(1-b)/b+6)/beta((a*(1-b))/b,(1-a)(1-b)/b)+
                                        choose(10,5)*beta((a(1-b))/b+5,(1-a)(1-b)/b+5)/beta((a*(1-b))/b,(1-a)(1-b)/b)+
                                        choose(10,6)*beta((a(1-b))/b+6,(1-a)(1-b)/b+4)/beta((a*(1-b))/b,(1-a)*(1-b)/b))
Deriv(F6,x)
#> Error in simplifications[[sym.name]]: wrong arguments for subsetting an environment

r[[2,1]]=c(F2)
r[[2,2]]=c(F6)
#> Error in r[[2, 2]] = c(F6): [[ ]] improper number of subscripts
r
#> [[1]]
#> [1] 0
#> 
#> [[2]]
#> [[2]][[1]]
#> function (a, b) 
#> (choose(10, 0) * beta((a(1 - b))/b + 0, (1 - a)(1 - b)/b + 10)/beta((a(1 - 
#>     b))/b, (1 - a)(1 - b)/b) + choose(10, 1) * beta((a(1 - b))/b + 
#>     1, (1 - a)(1 - b)/b + 9)/beta((a * (1 - b))/b, (1 - a)(1 - 
#>     b)/b) + choose(10, 2) * beta((a(1 - b))/b + 2, (1 - a)(1 - 
#>     b)/b + 8)/beta((a * (1 - b))/b, (1 - a)(1 - b)/b))
#> <environment: 0x7ffde4d4ab20>
#> 
#> 
#> [[3]]
#> [1] 0
#> 
#> [[4]]
#> [1] 0

timtimalo · January 5, 2018, 6:54pm

I use *s to represent multiplication. i.e., choose is multiplied by beta and choose will give a number but beta is function in "a" and "b" (which I am interest to calculate partial derivative with respect to them).
I fixed the name of matrix to be "mat" instead of "r" (to store the results of partial derivatives in it) in the second code I sent it.
Thank you