Отправка пакета CRAN: «Ошибка: C использование стека слишком близко к пределу» - PullRequest
2 голосов
/ 17 января 2020

Правильно: это проблема, с которой я столкнулся при отправке пакета R в CRAN. Поэтому я

  • не могу контролировать размер стека (поскольку проблема возникла на одной из платформ CRAN)
  • Я не могу предоставить воспроизводимый пример (поскольку я не знаю точных конфигураций в CRAN)

Проблема

При попытке отправить пакет cSEM.DGP в CRAN автоматический предварительный тест c (для Debian x86_64-p c - linux -gnu; не для Windows!) завершился неудачей с ПРИМЕЧАНИЕ: C stack usage 7975520 is too close to the limit.

Я знаю, что это вызвано функцией с тремя аргументами, тело которой имеет длину около 800 строк. Тело функции состоит из сложений и умножений этих аргументов. Это функция varzeta6(), которую вы найдете здесь (от строки 647 и далее).

Как я могу решить эту проблему?

То, что я не могу сделать:

  • предоставить воспроизводимый пример (по крайней мере, я бы не знал, как)
  • изменить размер стека

Вещи, о которых я думаю:

  • попробуйте разбить функцию на более мелкие части. Но я не знаю, как лучше всего это сделать.
  • как-то прекомпилировать? функция (если честно, я просто догадываюсь), чтобы CRAN не жаловался?

Дайте мне знать ваши идеи!

Подробности / Фон

Причина, почему varzeta6()varzeta4() / varzeta5() и даже более того varzeta7()) настолько длинны и R-неэффективны, что по существу они вставляются с помощью mathematica (после максимально возможного упрощения кода mathematica и его адаптации). быть действительным R-кодом). Следовательно, код ни в коем случае не оптимизирован для R (на что @MauritsEvers правильно указал).

Зачем нам математика? Потому что нам нужен общий вид корреляционной матрицы конструктивно подразумеваемой конструкции модели рекурсивного структурного уравнения, содержащей до 8 конструкций в зависимости от параметров уравнений модели. Кроме того, есть ограничения. Чтобы понять проблему, давайте возьмем систему из двух уравнений, которые могут быть рекурсивно решены:

  • Y2 = beta1 * Y1 + zeta1
  • Y3 = beta2 * Y1 + beta3 * Y2 + zeta2

Нам интересны ковариации: E (Y1 * Y2), E (Y1 * Y3) и E (Y2 * Y3) как функция от beta1, beta2 , бета3 при условии, что

  • E (Y1) = E (Y2) = E (Y3) = 0,
  • E (Y1 ^ 2) = E (Y2 ^ 2 ) = E (Y3 ^ 3) = 1
  • E (Yi * zeta_j) = 0 (при i = 1, 2, 3 и j = 1, 2)

Для такая простая модель, это довольно тривиально:

  • E (Y1 * Y2) = E (Y1 * (бета1 * Y1 + zeta1) = бета1 * E (Y1 ^ 2) + E (Y1 * zeta1) = beta1
  • E (Y1 * Y3) = E (Y1 * (beta2 * Y1 + beta3 * (beta1 * Y1 + zeta1) + zeta2) = beta2 + beta3 * beta1
  • E (Y2 * Y3) = ...

Но вы видите, как быстро это становится грязным, когда вы добавляете Y4, Y5, до Y8. В общем случае можно записать матрицу корреляции конструктивно подразумеваемой конструкции как (выражение на самом деле выглядит более сложным, потому что использование мы также допускаем до 5 экзогенных конструкций. Вот почему varzeta1() уже выглядит сложным. Но пока проигнорируйте это.):

  • V (Y) = (I - B) ^ - 1 В (дзета) (I - B) '^ - 1

где I - единичная матрица, а B - нижняя три angular матрица параметров модели (бета). V (дзета) является диагональной матрицей. Функции varzeta1(), varzeta2(), ..., varzeta7() вычисляют основные диагональные элементы. Поскольку мы ограничиваем Var (Yi) всегда равным 1, дисперсии Зетов следуют. Возьмем для примера уравнение Var (Y2) = beta1 ^ 2 * Var (Y1) + Var (zeta1) -> Var (zeta1) = 1 - beta1 ^ 2. Здесь это выглядит просто, но становится чрезвычайно сложным, когда мы берем дисперсию, скажем, шестого уравнения в такой цепочке рекурсивных уравнений, потому что Var (zeta6) зависит от всех предыдущих ковариаций между Y1, ..., Y5, которые сами по себе зависит от их соответствующих предыдущих ковариаций.

Хорошо, я не знаю, делает ли это что-то более ясным. Вот основной смысл:

  1. Код для varzeta1(), ..., varzeta7() является копией, вставленной из mathematica и, следовательно, не оптимизированной для R.
  2. Mathematica требуется потому что, насколько я знаю, R не может обрабатывать символьные c вычисления.
  3. Я мог бы оптимизировать R "вручную" (что крайне утомительно)
  4. Я думаю, что структура varzetaX() должен быть принят как дано. Поэтому возникает вопрос: могу ли я так или иначе использовать эту функцию?

Ответы [ 2 ]

3 голосов
/ 24 января 2020

Один из возможных подходов - попытаться убедить сопровождающих CRAN в том, что у вас нет простого способа решить проблему. Это NOTE, а не WARNING; Политика CRAN-репозитория гласит:

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

Итак, Вы могли бы рискнуть, что ваше аргументированное объяснение (в поле для комментариев на форме представления) убедит сопровождающих CRAN. В долгосрочной перспективе было бы лучше найти способ упростить вычисления, но, возможно, нет необходимости делать это перед отправкой в ​​CRAN.

1 голос
/ 24 января 2020

Это слишком длинный комментарий, но, надеюсь, это даст вам некоторые идеи по оптимизации кода для функций varzeta*; или, по крайней мере, это может дать вам пищу для размышлений.

Меня смущает несколько вещей:

  1. Все функции varzeta* имеют аргументы beta, gamma и phi, которые кажутся матрицами. Однако в varzeta1 вы не используете beta, но beta является первым аргументом функции.
  2. Я изо всех сил пытаюсь связать детали, которые вы даете в нижней части вашего поста, с кодом для varzeta* функций. Вы не объясняете, откуда взялись матрицы gamma и phi, и что они обозначают. Кроме того, учитывая, что beta являются оценками параметров модели, я не понимаю, почему beta должна быть матрицей.

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

Например, вы можете использовать crossprod и tcrossprod для вычисления перекрестных произведений, и %*% реализует умножение матриц.

Во-вторых, многие математические операции в R векторизованы. Я уже упоминал, что вы можете упростить

1 - gamma[1,1]^2 - gamma[1,2]^2 - gamma[1,3]^2 - gamma[1,4]^2 - gamma[1,5]^2

как

1 - sum(gamma[1, ]^2)

, так как оператор ^ векторизован.


Возможно, более фундаментально, это кажется что-то вроде проблемы XY, где это может помочь сделать шаг назад. Не зная полных деталей того, что вы пытаетесь смоделировать (как я уже сказал, я не могу связать данные, которые вы даете с кодом cSEM.DGP), я бы начал с изучения того, как решить рекурсивную SEM в R. на самом деле не вижу необходимости в Mathematica здесь. Как я уже говорил ранее, матричные операции очень стандартны в R; аналитическое решение набора рекурсивных уравнений также возможно в R. Так как вы, кажется, пришли из области Mathematica, было бы неплохо обсудить это с местным экспертом по R-кодированию.

Если вы должны использовать эти страшные varzeta* функции (и я действительно в этом сомневаюсь), можно переписать их в C ++, а затем скомпилировать с помощью Rcpp, чтобы превратить их в функции R. Возможно, это позволит избежать ограничения использования стека C?

...