Обсуждение статьи "Статистические распределения в MQL5 - берем лучшее из R и делаем быстрее" - страница 7

 
ivanivan_11:

я хоть и говорю о R, однако мой скилл вери смолл)) кто-то может проверить правильность кода?

если код верен, то можете проверить бенчмарк?

Код ошибочен.

Вы же замерили время компиляции функции, а не ее исполнение:

The function cmpfun compiles the body of a closure and returns a new closure with the same formals and the body replaced by the compiled body expression.


Доказательство ошибочности:

> library(microbenchmark)
> library(compiler)
> n <- 2000
> k <- seq(0,n,by=20)
> qqq<- function(xx) { a <- dbinom(k, n, pi/10, log = TRUE) }
> res <- microbenchmark(cmpfun(qqq))
> a
Ошибка: объект 'a' не найден

Если бы  во время бенчмарка запускалась функция qqq, то объект a получил бы расчитанные данные. Но вместо этого оказалось, что объект даже не был создан.

В результате бенчмарк посчитал время компиляции, а не время исполнения. В моем коде все верно - бенчмарк считает фактическое время исполнения и объект a создается с правильными данными.


И да, компиляция - достаточно затратный процесс: он в миллисекундах показан вместо микросекунд

> print(res)
Unit: milliseconds
        expr      min       lq     mean   median       uq      max neval
 cmpfun(qqq) 1.741577 1.783867 1.856366 1.812613 1.853096 2.697548   100
 

Ну и отдельный прикол, как вы в своем примере переопределили системную функцию q() - quit.

Выйти из R нельзя было никак :)

 
Dr.Trader:

Я к тому что mql компилятор во время компиляции уже знает все входные параметры. Ему достаточно во время компиляции всё посчитать, а при вызове скрипта - он просто возвращает уже заранее посчитанный результат. Видел какие-то статьи на хабре где сравнивали компиляторы на c++, там судя по анализу ассемблерного кода именно так происходит.

Да, может и активно этим пользуется. Вот примеры: https://www.mql5.com/ru/forum/58241

Но в данном случае это не пройдет - считать надо в полный рост из-за сложности, цикла и заполнения массива.

 
ivanivan_11:

если код верен, то можете проверить бенчмарк?

Нужно res <- microbenchmark(cmpfun(q)) заменить на res <- microbenchmark(q()). Но уже ранее скомпилированные библиотеки в байткод не перекомпилятся, у меня резульаты остались теже.

 

Renat Fatkhullin:
qqq<- function(xx) { a <- dbinom(k, n, pi/10, log = TRUE) }

"a" в таком случае будет локальной переменной, недоступной за пределами самой функции в любом случае. Но можно сделать так -
a <<- dbinom(k, n, pi/10, log = TRUE)
тогда она будет таки глобальной

 

Renat Fatkhullin:

Но в данном случае это не пройдет - считать надо в полный рост из-за сложности, цикла и заполнения массива.

понятно, скорость выполнения отличная тогда

 

Кстати, интерпретировать примитивный вызов a <- dbinom(k, n, pi/10, log = TRUE) с прямым проваливанием в ядро R с нативным исполнением (dbinom находится в r.dll) практически не стоит ничего.

Так что попытки скомпилировать этот вызов заведомо бесполезны.

 

Так как я много раз писал о быстродействии R, то позвольте вставить свои пять копеек.

 

Уважаемый Ренат!

Ваш пример вообще не о чем!

Вы взяли две аналогичные функции и делаете вывод вообще о быстродействии R.

Функции, которые Вы привели, вообще не отображают мощь и многообразие R.

Сравнивать надо вычислительно емкие операции.

Например, перемножение матриц...

Давайте померяем выражение на R

с <- a*b,

где a и b матрицы размером хотя бы 100*100. При этом Вы в своем коде проследите, чтобы R использовал MKL Intel. А это достигается просто установкой соответствующей версии R.

 

 

Если мы посмотрим на R, то горы кода, содержащего вычислительно емкие операции. Для их исполнения используются библиотеки, которые на данный момент максимальную эффективность.

И полезность Rв трейдинге не в функциях,которые Вы переписали (хотя они также необходимы), а в моделях. В одном из ответов Вам я упоминал пакет caret. Посмотрите, что это такое.... Вот реализация какой-либо практически полезной торговой модели  в рамках этого пакета и на мкл и даст ответ

 

Кроме этого следует не забывать, что загрузка всех ядер компа - это рутина для R. Кроме этого можно загрузить соседние компы в локальной сети.  

 

ПС.

Для меня является сомнительной идея сравнения быстродействия МКЛ и R: эти две системы имеют совершенно разные предметные области 

 

СанСаныч, все проверим и выпустим бенчмарк. Но сначала допишем функционал.

Тест был обоснован и он сразу вскрыл проблему. Теоретическое обоснование я представил и уверен, что системный оверхед у R сохранится по подавляющему объему функционала.

Матрицы и мы перемножить можем так, что и Интел проиграет. Это давно уже не rocket science, а Интел(вернее, такие же сторонние программисты в рамках принадлежности к компании) не чемпион в мифическом знании своих процессоров.

 
СанСаныч Фоменко:

Так как я много раз писал о быстродействии R, то позвольте вставить свои пять копеек.

 .................

Сан-Санычу и другим ребятам.

Сан-Саныч, Вы же знаете как я Вас уважаю ... ((С) Катаев и Файнзильберг, известные как "Ильф и Петров"), несмотря на некоторые Ваши пост-советские приколы тут.

Позвольте я Вам поясню кое-что важное:

1). Основная работа программиста - это не писать программы, а ЧИТАТЬ программы, в частности свою. Любой программист 95...99% своего времени сидит и пялится в монитор. Разве он ПИШЕТ программу? Нет, он её в основном ЧИТАЕТ. Поэтому чем ближе к естественному языку будет то, что он читает на экране, то есть к тому, чему его учили мама, папа, бабушка, учительница в школе, - тем эффективнее он будет расшифровывать эти языко-пободные кракозябры на экране и втыкать соответствие между алгоритмом и своей программой.

2). Для целей пункта (1) нет ничего лучше в среднем, чем язык Си.  именно поэтому, например мне лично (а также 2-3 ответственных, и не очень, человека) удалось написать проект объёмом в 700+ подпрограмм, на Си, MQL4, CUDA.... И всё работает.

3). С точки зрения пункта (1) объектно-ориентированный вариант Си, то есть Си++ намного хуже. (Но об этом в другой раз).

4). Полная совместимость классического Си и MQL4 просто бесценна. Перенос процедуры туда-обратно занимает полминуты.

5). Основное достоинство Си+MQL4 - это CLARITY. То есть понятность и прозрачность всего, что есть на экране у программиста.

Если сравнить Си-MQL4 с этой Вашей R, то здесь надо смотреть не на скорость и объём понаписанного кода, а на CLARITY текста. То есть на его понятность. Иначе программист будет сутками пялиться на экран, в тщетных попытках понять, что делает прога, какие у неё там параметры, почему автор их так назвал,  и вообще почему программер сделал именно так, а не иначе. Тут не скорость программы важна, а правильность её работы, и скорость ЕЁ ПРИМЕНИМОСТИ у конечного программиста.

С этой точки зрения то, что сделали Метаквоты - конечно огромная поддержка тем, кто захочет вставить статистику в свои советники. Тут даже нечего сравнивать по степени простоты и понятности функций. А это важно. Особенно если у Вас тонкие расчёты (а Форекс и трейдинг вообще требует тонких расчётов).

Давайте сравним.

Вот как выглядит на Си - MQL4 функция интегрирования:

//__________________________________________________________________
//                                                                                                                                                      |
//               Integral_Simpson_3_Points_Lite ()                                      |
//_______________________________________________________|
#define FACTOR_1_3      (1.0 / 3.0)

double Integral_Simpson_3_Points_Lite
(
   double       & Price_Array[],
   int          Period_Param,
   int          Last_Bar
)
{
        double  Sum, Sum_Full_Cycles, Sum_Tail, Sum_Total;
        int             i, j;
        int             Quant;
        int             Full_Cycles;
        int             Tail_Limit;
        int             Tail_Start;

        //..................................................................
        if (Last_Bar < 0)
        {
                Last_Bar = 0;
        }
        //.........................................
        if (Period_Param <= 1)
        {
                return (0.0);
        }
        //.................................................................
        if (Period_Param == 2)
        {
                return (0.5 * (Price_Array[Last_Bar] + Price_Array[Last_Bar + 1]));
        }
        //=============================================================================================
        //=============================================================================================
        Quant = 3;
        Full_Cycles = (Period_Param - 1) / (Quant - 1);
        Tail_Start = Full_Cycles * (Quant - 1);
        Tail_Limit = Period_Param - Tail_Start;
        //...........................................................
        j = Last_Bar;

        Sum = 0.0;
        for (i = 0; i < Full_Cycles; i ++)
        {
                //.........................................................................
                Sum += Price_Array[j];
                Sum += 4.0 * Price_Array[j + 1];
                Sum += Price_Array[j + 2];
                j = j + (Quant - 1);
        }
        //...........................................................
        Sum_Full_Cycles = Sum * FACTOR_1_3;
        Sum_Tail = Integral_Trapezoid_Lite (Price_Array,
                                            Tail_Limit,
                                            Last_Bar + Tail_Start);
        //.........................................................................
        Sum_Total = Sum_Full_Cycles + Sum_Tail;
        //...............................................................
        return (Sum_Total) ;

}
 

Буду писать частями, так легче писать.

Там есть внутри функция интегрирования трапецией:

//__________________________________________________________________
//                                                                                                                                                      |
//               Integral_Trapezoid_Lite ()                                                     |
//_______________________________________________________|
double Integral_Trapezoid_Lite
(
   double       & Price_Array[],
   int          Period_Param,
   int          Last_Bar
)
{
        double  Sum;
        int             i;
        int             Price_Index ;
        //.........................................
        if (Last_Bar < 0)
        {
                Last_Bar = 0;
        }
        //.........................................
        if (Period_Param <= 1)
        {
                return (0.0);
        }
        //.................................................................
        if (Period_Param == 2)
        {
                return (0.5 * (Price_Array[Last_Bar] + Price_Array[Last_Bar + 1]));
        }
        //..................................................................
        //..................................................................
        Sum = 0.0;
        for (i = 0; i < Period_Param; i++)
        {
                Price_Index = Last_Bar + i;
                if (Price_Index < 0)
                {
                        break;
                }
                //........................................
                if ((i == 0) || (i == (Period_Param - 1)))
                {
                        Sum = Sum + Price_Array[i] * 0.5;
                }
                else
                {
                        Sum = Sum + Price_Array[i];
                }
        }
        //...............................................................
        return (Sum) ;
}
 

Всё абсолютно ясно и понятно. И что важно, это работает всегда и работает хорошо, то есть с низкой погрешностью даже в MT4-MQL4, что экономит кучу времени.

Но если Вы захотите узнать почему у Вас вылазят непонятные погрешности при работе в R, или если Вы просто захотите разобраться какие там параметры в процедуре интегрирования или какой именно метод интегрирования они там запрограммировали, то Вы увидите следующее (прости Господи меня за выкладывание такого - незрелым отрокам программирования):

http://www.netlib.org/quadpack/

Это только заголовок функции, изначально написанной на Фортране. Основной текст будет позже. Это оригинал программы, которая используется в пакете R для интегрирования.

Что здесь можно понять, скажите мне?

      subroutine qagse(f,a,b,epsabs,epsrel,limit,result,abserr,neval,
     *   ier,alist,blist,rlist,elist,iord,last)
c***begin prologue  qagse
c***date written   800101   (yymmdd)
c***revision date  830518   (yymmdd)
c***category no.  h2a1a1
c***keywords  automatic integrator, general-purpose,
c             (end point) singularities, extrapolation,
c             globally adaptive
c***author  piessens,robert,appl. math. & progr. div. - k.u.leuven
c           de doncker,elise,appl. math. & progr. div. - k.u.leuven
c***purpose  the routine calculates an approximation result to a given
c            definite integral i = integral of f over (a,b),
c            hopefully satisfying following claim for accuracy
c            abs(i-result).le.max(epsabs,epsrel*abs(i)).
c***description
c
c        computation of a definite integral
c        standard fortran subroutine
c        real version
c
c        parameters
c         on entry
c            f      - real
c                     function subprogram defining the integrand
c                     function f(x). the actual name for f needs to be
c                     declared e x t e r n a l in the driver program.
c
c            a      - real
c                     lower limit of integration
c
c            b      - real
c                     upper limit of integration
c
c            epsabs - real
c                     absolute accuracy requested
c            epsrel - real
c                     relative accuracy requested
c                     if  epsabs.le.0
c                     and epsrel.lt.max(50*rel.mach.acc.,0.5d-28),
c                     the routine will end with ier = 6.
c
c            limit  - integer
c                     gives an upperbound on the number of subintervals
c                     in the partition of (a,b)
c
c         on return
c            result - real
c                     approximation to the integral
c
c            abserr - real
c                     estimate of the modulus of the absolute error,
c                     which should equal or exceed abs(i-result)
c
c            neval  - integer
c                     number of integrand evaluations
c
c            ier    - integer
c                     ier = 0 normal and reliable termination of the
c                             routine. it is assumed that the requested
c                             accuracy has been achieved.
c                     ier.gt.0 abnormal termination of the routine
c                             the estimates for integral and error are
c                             less reliable. it is assumed that the
c                             requested accuracy has not been achieved.
c            error messages
c                         = 1 maximum number of subdivisions allowed
c                             has been achieved. one can allow more sub-
c                             divisions by increasing the value of limit
c                             (and taking the according dimension
c                             adjustments into account). however, if
c                             this yields no improvement it is advised
c                             to analyze the integrand in order to
c                             determine the integration difficulties. if
c                             the position of a local difficulty can be
c                             determined (e.g. singularity,
c                             discontinuity within the interval) one
c                             will probably gain from splitting up the
c                             interval at this point and calling the
c                             integrator on the subranges. if possible,
c                             an appropriate special-purpose integrator
c                             should be used, which is designed for
c                             handling the type of difficulty involved.
c                         = 2 the occurrence of roundoff error is detec-
c                             ted, which prevents the requested
c                             tolerance from being achieved.
c                             the error may be under-estimated.
c                         = 3 extremely bad integrand behaviour
c                             occurs at some points of the integration
c                             interval.
c                         = 4 the algorithm does not converge.
c                             roundoff error is detected in the
c                             extrapolation table.
c                             it is presumed that the requested
c                             tolerance cannot be achieved, and that the
c                             returned result is the best which can be
c                             obtained.
c                         = 5 the integral is probably divergent, or
c                             slowly convergent. it must be noted that
c                             divergence can occur with any other value
c                             of ier.
c                         = 6 the input is invalid, because
c                             epsabs.le.0 and
c                             epsrel.lt.max(50*rel.mach.acc.,0.5d-28).
c                             result, abserr, neval, last, rlist(1),
c                             iord(1) and elist(1) are set to zero.
c                             alist(1) and blist(1) are set to a and b
c                             respectively.
c
c            alist  - real
c                     vector of dimension at least limit, the first
c                      last  elements of which are the left end points
c                     of the subintervals in the partition of the
c                     given integration range (a,b)
c
c            blist  - real
c                     vector of dimension at least limit, the first
c                      last  elements of which are the right end points
c                     of the subintervals in the partition of the given
c                     integration range (a,b)
c
c            rlist  - real
c                     vector of dimension at least limit, the first
c                      last  elements of which are the integral
c                     approximations on the subintervals
c
c            elist  - real
c                     vector of dimension at least limit, the first
c                      last  elements of which are the moduli of the
c                     absolute error estimates on the subintervals
c
c            iord   - integer
c                     vector of dimension at least limit, the first k
c                     elements of which are pointers to the
c                     error estimates over the subintervals,
c                     such that elist(iord(1)), ..., elist(iord(k))
c                     form a decreasing sequence, with k = last
c                     if last.le.(limit/2+2), and k = limit+1-last
c                     otherwise
c
c            last   - integer
c                     number of subintervals actually produced in the
c                     subdivision process
c
c***references  (none)
c***routines called  qelg,qk21,qpsrt,r1mach
c***end prologue  qagse
c
      real a,abseps,abserr,alist,area,area1,area12,area2,a1,
     *  a2,b,blist,b1,b2,correc,defabs,defab1,defab2,r1mach,
     *  dres,elist,epmach,epsabs,epsrel,erlarg,erlast,errbnd,
     *  errmax,error1,error2,erro12,errsum,ertest,f,oflow,resabs,
     *  reseps,result,res3la,rlist,rlist2,small,uflow
      integer id,ier,ierro,iord,iroff1,iroff2,iroff3,jupbnd,k,ksgn,
     *  ktmin,last,limit,maxerr,neval,nres,nrmax,numrl2
      logical extrap,noext
c
      dimension alist(limit),blist(limit),elist(limit),iord(limit),
     * res3la(3),rlist(limit),rlist2(52)
c
      external f
c
c            the dimension of rlist2 is determined by the value of
c            limexp in subroutine qelg (rlist2 should be of dimension
c            (limexp+2) at least).
c
c            list of major variables
c            -----------------------
c
c           alist     - list of left end points of all subintervals
c                       considered up to now
c           blist     - list of right end points of all subintervals
c                       considered up to now
c           rlist(i)  - approximation to the integral over
c                       (alist(i),blist(i))
c           rlist2    - array of dimension at least limexp+2
c                       containing the part of the epsilon table
c                       which is still needed for further computations
c           elist(i)  - error estimate applying to rlist(i)
c           maxerr    - pointer to the interval with largest error
c                       estimate
c           errmax    - elist(maxerr)
c           erlast    - error on the interval currently subdivided
c                       (before that subdivision has taken place)
c           area      - sum of the integrals over the subintervals
c           errsum    - sum of the errors over the subintervals
c           errbnd    - requested accuracy max(epsabs,epsrel*
c                       abs(result))
c           *****1    - variable for the left interval
c           *****2    - variable for the right interval
c           last      - index for subdivision
c           nres      - number of calls to the extrapolation routine
c           numrl2    - number of elements currently in rlist2. if an
c                       appropriate approximation to the compounded
c                       integral has been obtained it is put in
c                       rlist2(numrl2) after numrl2 has been increased
c                       by one.
c           small     - length of the smallest interval considered
c                       up to now, multiplied by 1.5
c           erlarg    - sum of the errors over the intervals larger
c                       than the smallest interval considered up to now
c           extrap    - logical variable denoting that the routine
c                       is attempting to perform extrapolation
c                       i.e. before subdividing the smallest interval
c                       we try to decrease the value of erlarg.
c           noext     - logical variable denoting that extrapolation
c                       is no longer allowed (true value)
c
c            machine dependent constants
c            ---------------------------
c
c           epmach is the largest relative spacing.
c           uflow is the smallest positive magnitude.
c           oflow is the largest positive magnitude.
c
c***first executable statement  qagse
      epmach = r1mach(4)
c
c            test on validity of parameters
c            ------------------------------