[R] multiply of two expressions

Gabor Grothendieck ggrothendieck at gmail.com
Tue May 6 04:56:22 CEST 2014


On Mon, May 5, 2014 at 10:16 PM, David Winsemius <dwinsemius at comcast.net> wrote:
>
> On May 5, 2014, at 6:05 PM, Gabor Grothendieck wrote:
>
>> On Mon, May 5, 2014 at 9:43 AM, Niloofar.Javanrouh
>> <javanrouh_n at yahoo.com> wrote:
>>>
>>>
>>> hello,
>>> i want to differentiate of L with respect to b
>>> when:
>>>
>>> L= k*ln (k/(k+mu)) + sum(y) * ln (1-(k/mu+k))   #(negative binomial ln likelihood)
>>> and
>>> ln(mu/(mu+k)) = a+bx   #link function
>>>
>>> how can i do it in R?
>>
>> Try this.  First we solve for 'mu' in terms of the other variables
>> using the link equation:
>>
>>> library(Ryacas)
>>> k <- Sym("k")
>>> mu <- Sym("mu")
>>> y <- Sym("y")
>>> L <- Sym("L")
>>> a <- Sym("a")
>>> b <- Sym("b")
>>> x <- Sym("x")
>>> sumy <- Sym("sumy")
>>> Solve(log(mu/(mu+k)) == a+b*x, "mu")
>> expression(list(mu == k * exp(a + b * x)/(1 - exp(a + b * x))))
>>>
>>
>> Now in 'k*log(k/(k+mu)) + sumy * log(1-(k/mu+k))' substitute 'k *
>> exp(a + b * x)/(1 - exp(a + b * x)))' for 'mu' using 'Subst' and take
>> the derivative with respect to 'a' using 'deriv':
>>
>>> s <- Subst(k*log(k/(k+mu)) + sumy * log(1-(k/mu+k)), mu,
>> +       k * exp(a + b * x)/(1 - exp(a + b * x)))
>>> deriv(s, a)
>> expression(sumy * (k * ((1 - exp(a + b * x)) * (k * exp(a + b *
>>    x)) + k * exp(a + b * x)^2))/((1 - exp(a + b * x))^2 * (k *
>>    exp(a + b * x)/(1 - exp(a + b * x)))^2 * (1 - (k * (1 - exp(a +
>>    b * x))/(k * exp(a + b * x)) + k))) - k * ((k + k * exp(a +
>>    b * x)/(1 - exp(a + b * x))) * (k * ((1 - exp(a + b * x)) *
>>    (k * exp(a + b * x)) + k * exp(a + b * x)^2)))/((1 - exp(a +
>>    b * x))^2 * ((k + k * exp(a + b * x)/(1 - exp(a + b * x)))^2 *
>>    k)))
>>
>
> But can we maybe get the Taylor series approximation to first or second order?
>



> Taylor(d, a, 0, 2)
expression(sumy * (k * ((1 - exp(b * x)) * (k * exp(b * x)) +
    k * exp(b * x)^2))/((1 - exp(b * x))^2 * (k * exp(b * x)/(1 -
    exp(b * x)))^2 * (1 - (k * (1 - exp(b * x))/(k * exp(b *
    x)) + k))) - k * ((k + k * exp(b * x)/(1 - exp(b * x))) *
    (k * ((1 - exp(b * x)) * (k * exp(b * x)) + k * exp(b * x)^2)))/((1 -
    exp(b * x))^2 * ((k + k * exp(b * x)/(1 - exp(b * x)))^2 *
    k)) + (((1 - exp(b * x))^2 * (k * exp(b * x)/(1 - exp(b *
    x)))^2 * (1 - (k * (1 - exp(b * x))/(k * exp(b * x)) + k)) *
    (sumy * (k * ((1 - exp(b * x)) * (k * exp(b * x)) - k * exp(b *
        x)^2 + k * (2 * exp(b * x)^2)))) - sumy * (k * ((1 -
    exp(b * x)) * (k * exp(b * x)) + k * exp(b * x)^2)) * ((1 -
    exp(b * x))^2 * (k * exp(b * x)/(1 - exp(b * x)))^2 * (k *
    ((1 - exp(b * x)) * (k * exp(b * x)) + k * exp(b * x)^2))/((1 -
    exp(b * x))^2 * (k * exp(b * x)/(1 - exp(b * x)))^2) + ((1 -
    exp(b * x))^2 * (k * exp(b * x) * (2 * ((1 - exp(b * x)) *
    (k * exp(b * x)) + k * exp(b * x)^2)))/(1 - exp(b * x))^3 +
    -2 * exp(b * x) * (1 - exp(b * x)) * (k * exp(b * x)/(1 -
        exp(b * x)))^2) * (1 - (k * (1 - exp(b * x))/(k * exp(b *
    x)) + k))))/((1 - exp(b * x))^2 * (k * exp(b * x)/(1 - exp(b *
    x)))^2 * (1 - (k * (1 - exp(b * x))/(k * exp(b * x)) + k)))^2 -
    ((1 - exp(b * x))^2 * ((k + k * exp(b * x)/(1 - exp(b * x)))^2 *
        k) * (k * ((k + k * exp(b * x)/(1 - exp(b * x))) * (k *
        ((1 - exp(b * x)) * (k * exp(b * x)) - k * exp(b * x)^2 +
            k * (2 * exp(b * x)^2))) + k * ((1 - exp(b * x)) *
        (k * exp(b * x)) + k * exp(b * x)^2)^2/(1 - exp(b * x))^2)) -
        k * ((k + k * exp(b * x)/(1 - exp(b * x))) * (k * ((1 -
            exp(b * x)) * (k * exp(b * x)) + k * exp(b * x)^2))) *
            ((1 - exp(b * x))^2 * (k * ((k + k * exp(b * x)/(1 -
                exp(b * x))) * (2 * ((1 - exp(b * x)) * (k *
                exp(b * x)) + k * exp(b * x)^2))))/(1 - exp(b *
                x))^2 + -2 * exp(b * x) * (1 - exp(b * x)) *
                ((k + k * exp(b * x)/(1 - exp(b * x)))^2 * k)))/((1 -
        exp(b * x))^2 * ((k + k * exp(b * x)/(1 - exp(b * x)))^2 *
        k))^2) * a + a^2 * ((((1 - exp(b * x))^2 * (k * exp(b *
    x)/(1 - exp(b * x)))^2 * (1 - (k * (1 - exp(b * x))/(k *
    exp(b * x)) + k)))^2 * ((1 - exp(b * x))^2 * (k * exp(b *
    x)/(1 - exp(b * x)))^2 * (1 - (k * (1 - exp(b * x))/(k *
    exp(b * x)) + k)) * (sumy * (k * ((1 - exp(b * x)) * (k *
    exp(b * x)) - k * exp(b * x)^2 - k * (2 * exp(b * x)^2) +
    k * (4 * exp(b * x)^2)))) + ((1 - exp(b * x))^2 * (k * exp(b *
    x)/(1 - exp(b * x)))^2 * (k * ((1 - exp(b * x)) * (k * exp(b *
    x)) + k * exp(b * x)^2))/((1 - exp(b * x))^2 * (k * exp(b *
    x)/(1 - exp(b * x)))^2) + ((1 - exp(b * x))^2 * (k * exp(b *
    x) * (2 * ((1 - exp(b * x)) * (k * exp(b * x)) + k * exp(b *
    x)^2)))/(1 - exp(b * x))^3 + -2 * exp(b * x) * (1 - exp(b *
    x)) * (k * exp(b * x)/(1 - exp(b * x)))^2) * (1 - (k * (1 -
    exp(b * x))/(k * exp(b * x)) + k))) * (sumy * (k * ((1 -
    exp(b * x)) * (k * exp(b * x)) - k * exp(b * x)^2 + k * (2 *
    exp(b * x)^2)))) - (sumy * (k * ((1 - exp(b * x)) * (k *
    exp(b * x)) + k * exp(b * x)^2)) * (((1 - exp(b * x))^2 *
    (k * exp(b * x)/(1 - exp(b * x)))^2 * ((1 - exp(b * x))^2 *
    (k * exp(b * x)/(1 - exp(b * x)))^2 * (k * ((1 - exp(b *
    x)) * (k * exp(b * x)) - k * exp(b * x)^2 + k * (2 * exp(b *
    x)^2))) + ((1 - exp(b * x))^2 * (k * exp(b * x) * (2 * ((1 -
    exp(b * x)) * (k * exp(b * x)) + k * exp(b * x)^2)))/(1 -
    exp(b * x))^3 + -2 * exp(b * x) * (1 - exp(b * x)) * (k *
    exp(b * x)/(1 - exp(b * x)))^2) * (k * ((1 - exp(b * x)) *
    (k * exp(b * x)) + k * exp(b * x)^2))) - (1 - exp(b * x))^2 *
    (k * exp(b * x)/(1 - exp(b * x)))^2 * (k * ((1 - exp(b *
    x)) * (k * exp(b * x)) + k * exp(b * x)^2)) * ((1 - exp(b *
    x))^2 * (k * exp(b * x) * (2 * ((1 - exp(b * x)) * (k * exp(b *
    x)) + k * exp(b * x)^2)))/(1 - exp(b * x))^3 + -2 * exp(b *
    x) * (1 - exp(b * x)) * (k * exp(b * x)/(1 - exp(b * x)))^2))/((1 -
    exp(b * x))^2 * (k * exp(b * x)/(1 - exp(b * x)))^2)^2 +
    (((1 - exp(b * x))^2 * (k * exp(b * x) * (2 * ((1 - exp(b *
        x)) * (k * exp(b * x)) + k * exp(b * x)^2)))/(1 - exp(b *
        x))^3 + -2 * exp(b * x) * (1 - exp(b * x)) * (k * exp(b *
        x)/(1 - exp(b * x)))^2) * (k * ((1 - exp(b * x)) * (k *
        exp(b * x)) + k * exp(b * x)^2))/((1 - exp(b * x))^2 *
        (k * exp(b * x)/(1 - exp(b * x)))^2) + (((1 - exp(b *
        x))^3 * ((1 - exp(b * x))^2 * (k * exp(b * x) * (2 *
        ((1 - exp(b * x)) * (k * exp(b * x)) - k * exp(b * x)^2 +
            k * (2 * exp(b * x)^2))) + k * exp(b * x) * (2 *
        ((1 - exp(b * x)) * (k * exp(b * x)) + k * exp(b * x)^2))) +
        -2 * exp(b * x) * (1 - exp(b * x)) * (k * exp(b * x) *
            (2 * ((1 - exp(b * x)) * (k * exp(b * x)) + k * exp(b *
                x)^2)))) - (1 - exp(b * x))^2 * (k * exp(b *
        x) * (2 * ((1 - exp(b * x)) * (k * exp(b * x)) + k *
        exp(b * x)^2))) * (-3 * exp(b * x) * (1 - exp(b * x))^2))/(1 -
        exp(b * x))^6 + (-2 * exp(b * x) * (1 - exp(b * x)) *
        (k * exp(b * x) * (2 * ((1 - exp(b * x)) * (k * exp(b *
            x)) + k * exp(b * x)^2)))/(1 - exp(b * x))^3 + (-2 *
        exp(b * x) * (1 - exp(b * x)) - -2 * exp(b * x)^2) *
        (k * exp(b * x)/(1 - exp(b * x)))^2)) * (1 - (k * (1 -
        exp(b * x))/(k * exp(b * x)) + k)))) + sumy * (k * ((1 -
    exp(b * x)) * (k * exp(b * x)) - k * exp(b * x)^2 + k * (2 *
    exp(b * x)^2))) * ((1 - exp(b * x))^2 * (k * exp(b * x)/(1 -
    exp(b * x)))^2 * (k * ((1 - exp(b * x)) * (k * exp(b * x)) +
    k * exp(b * x)^2))/((1 - exp(b * x))^2 * (k * exp(b * x)/(1 -
    exp(b * x)))^2) + ((1 - exp(b * x))^2 * (k * exp(b * x) *
    (2 * ((1 - exp(b * x)) * (k * exp(b * x)) + k * exp(b * x)^2)))/(1 -
    exp(b * x))^3 + -2 * exp(b * x) * (1 - exp(b * x)) * (k *
    exp(b * x)/(1 - exp(b * x)))^2) * (1 - (k * (1 - exp(b *
    x))/(k * exp(b * x)) + k))))) - ((1 - exp(b * x))^2 * (k *
    exp(b * x)/(1 - exp(b * x)))^2 * (1 - (k * (1 - exp(b * x))/(k *
    exp(b * x)) + k)) * (sumy * (k * ((1 - exp(b * x)) * (k *
    exp(b * x)) - k * exp(b * x)^2 + k * (2 * exp(b * x)^2)))) -
    sumy * (k * ((1 - exp(b * x)) * (k * exp(b * x)) + k * exp(b *
        x)^2)) * ((1 - exp(b * x))^2 * (k * exp(b * x)/(1 - exp(b *
        x)))^2 * (k * ((1 - exp(b * x)) * (k * exp(b * x)) +
        k * exp(b * x)^2))/((1 - exp(b * x))^2 * (k * exp(b *
        x)/(1 - exp(b * x)))^2) + ((1 - exp(b * x))^2 * (k *
        exp(b * x) * (2 * ((1 - exp(b * x)) * (k * exp(b * x)) +
        k * exp(b * x)^2)))/(1 - exp(b * x))^3 + -2 * exp(b *
        x) * (1 - exp(b * x)) * (k * exp(b * x)/(1 - exp(b *
        x)))^2) * (1 - (k * (1 - exp(b * x))/(k * exp(b * x)) +
        k)))) * (2 * ((1 - exp(b * x))^2 * (k * exp(b * x)/(1 -
    exp(b * x)))^2 * (k * ((1 - exp(b * x)) * (k * exp(b * x)) +
    k * exp(b * x)^2))/((1 - exp(b * x))^2 * (k * exp(b * x)/(1 -
    exp(b * x)))^2) + ((1 - exp(b * x))^2 * (k * exp(b * x) *
    (2 * ((1 - exp(b * x)) * (k * exp(b * x)) + k * exp(b * x)^2)))/(1 -
    exp(b * x))^3 + -2 * exp(b * x) * (1 - exp(b * x)) * (k *
    exp(b * x)/(1 - exp(b * x)))^2) * (1 - (k * (1 - exp(b *
    x))/(k * exp(b * x)) + k))) * ((1 - exp(b * x))^2 * (k *
    exp(b * x)/(1 - exp(b * x)))^2 * (1 - (k * (1 - exp(b * x))/(k *
    exp(b * x)) + k)))))/((1 - exp(b * x))^2 * (k * exp(b * x)/(1 -
    exp(b * x)))^2 * (1 - (k * (1 - exp(b * x))/(k * exp(b *
    x)) + k)))^4 - (((1 - exp(b * x))^2 * ((k + k * exp(b * x)/(1 -
    exp(b * x)))^2 * k))^2 * ((1 - exp(b * x))^2 * ((k + k *
    exp(b * x)/(1 - exp(b * x)))^2 * k) * (k * ((k + k * exp(b *
    x)/(1 - exp(b * x))) * (k * ((1 - exp(b * x)) * (k * exp(b *
    x)) - k * exp(b * x)^2 - k * (2 * exp(b * x)^2) + k * (4 *
    exp(b * x)^2))) + k * ((1 - exp(b * x)) * (k * exp(b * x)) -
    k * exp(b * x)^2 + k * (2 * exp(b * x)^2)) * ((1 - exp(b *
    x)) * (k * exp(b * x)) + k * exp(b * x)^2)/(1 - exp(b * x))^2 +
    ((1 - exp(b * x))^2 * (k * (2 * ((1 - exp(b * x)) * (k *
        exp(b * x)) - k * exp(b * x)^2 + k * (2 * exp(b * x)^2)) *
        ((1 - exp(b * x)) * (k * exp(b * x)) + k * exp(b * x)^2))) -
        k * ((1 - exp(b * x)) * (k * exp(b * x)) + k * exp(b *
            x)^2)^2 * (-2 * exp(b * x) * (1 - exp(b * x))))/(1 -
        exp(b * x))^4)) + ((1 - exp(b * x))^2 * (k * ((k + k *
    exp(b * x)/(1 - exp(b * x))) * (2 * ((1 - exp(b * x)) * (k *
    exp(b * x)) + k * exp(b * x)^2))))/(1 - exp(b * x))^2 + -2 *
    exp(b * x) * (1 - exp(b * x)) * ((k + k * exp(b * x)/(1 -
    exp(b * x)))^2 * k)) * (k * ((k + k * exp(b * x)/(1 - exp(b *
    x))) * (k * ((1 - exp(b * x)) * (k * exp(b * x)) - k * exp(b *
    x)^2 + k * (2 * exp(b * x)^2))) + k * ((1 - exp(b * x)) *
    (k * exp(b * x)) + k * exp(b * x)^2)^2/(1 - exp(b * x))^2)) -
    (k * ((k + k * exp(b * x)/(1 - exp(b * x))) * (k * ((1 -
        exp(b * x)) * (k * exp(b * x)) + k * exp(b * x)^2))) *
        (((1 - exp(b * x))^2 * ((1 - exp(b * x))^2 * (k * ((k +
            k * exp(b * x)/(1 - exp(b * x))) * (2 * ((1 - exp(b *
            x)) * (k * exp(b * x)) - k * exp(b * x)^2 + k * (2 *
            exp(b * x)^2))) + 2 * ((1 - exp(b * x)) * (k * exp(b *
            x)) + k * exp(b * x)^2)^2/(1 - exp(b * x))^2)) +
            -2 * exp(b * x) * (1 - exp(b * x)) * (k * ((k + k *
                exp(b * x)/(1 - exp(b * x))) * (2 * ((1 - exp(b *
                x)) * (k * exp(b * x)) + k * exp(b * x)^2))))) -
            (1 - exp(b * x))^2 * (k * ((k + k * exp(b * x)/(1 -
                exp(b * x))) * (2 * ((1 - exp(b * x)) * (k *
                exp(b * x)) + k * exp(b * x)^2)))) * (-2 * exp(b *
                x) * (1 - exp(b * x))))/(1 - exp(b * x))^4 +
            (-2 * exp(b * x) * (1 - exp(b * x)) * (k * ((k +
                k * exp(b * x)/(1 - exp(b * x))) * (2 * ((1 -
                exp(b * x)) * (k * exp(b * x)) + k * exp(b *
                x)^2))))/(1 - exp(b * x))^2 + (-2 * exp(b * x) *
                (1 - exp(b * x)) - -2 * exp(b * x)^2) * ((k +
                k * exp(b * x)/(1 - exp(b * x)))^2 * k))) + k *
        ((k + k * exp(b * x)/(1 - exp(b * x))) * (k * ((1 - exp(b *
            x)) * (k * exp(b * x)) - k * exp(b * x)^2 + k * (2 *
            exp(b * x)^2))) + k * ((1 - exp(b * x)) * (k * exp(b *
            x)) + k * exp(b * x)^2)^2/(1 - exp(b * x))^2) * ((1 -
        exp(b * x))^2 * (k * ((k + k * exp(b * x)/(1 - exp(b *
        x))) * (2 * ((1 - exp(b * x)) * (k * exp(b * x)) + k *
        exp(b * x)^2))))/(1 - exp(b * x))^2 + -2 * exp(b * x) *
        (1 - exp(b * x)) * ((k + k * exp(b * x)/(1 - exp(b *
        x)))^2 * k)))) - ((1 - exp(b * x))^2 * ((k + k * exp(b *
    x)/(1 - exp(b * x)))^2 * k) * (k * ((k + k * exp(b * x)/(1 -
    exp(b * x))) * (k * ((1 - exp(b * x)) * (k * exp(b * x)) -
    k * exp(b * x)^2 + k * (2 * exp(b * x)^2))) + k * ((1 - exp(b *
    x)) * (k * exp(b * x)) + k * exp(b * x)^2)^2/(1 - exp(b *
    x))^2)) - k * ((k + k * exp(b * x)/(1 - exp(b * x))) * (k *
    ((1 - exp(b * x)) * (k * exp(b * x)) + k * exp(b * x)^2))) *
    ((1 - exp(b * x))^2 * (k * ((k + k * exp(b * x)/(1 - exp(b *
        x))) * (2 * ((1 - exp(b * x)) * (k * exp(b * x)) + k *
        exp(b * x)^2))))/(1 - exp(b * x))^2 + -2 * exp(b * x) *
        (1 - exp(b * x)) * ((k + k * exp(b * x)/(1 - exp(b *
        x)))^2 * k))) * (2 * ((1 - exp(b * x))^2 * (k * ((k +
    k * exp(b * x)/(1 - exp(b * x))) * (2 * ((1 - exp(b * x)) *
    (k * exp(b * x)) + k * exp(b * x)^2))))/(1 - exp(b * x))^2 +
    -2 * exp(b * x) * (1 - exp(b * x)) * ((k + k * exp(b * x)/(1 -
        exp(b * x)))^2 * k)) * ((1 - exp(b * x))^2 * ((k + k *
    exp(b * x)/(1 - exp(b * x)))^2 * k))))/((1 - exp(b * x))^2 *
    ((k + k * exp(b * x)/(1 - exp(b * x)))^2 * k))^4)/2)



More information about the R-help mailing list