pseudobot.saty

@require: stdjabook
@require: code
@require: itemize
@require: tabular
@require: math
@import: local

document (|
  title = {手動botのツイート間隔};
  author = {hsjoihs};
  show-title = true;
  show-toc = false;
|) '<
  +section{概要}<
    +p{
      およそ1時間おきにツイートする手動botを運用していたが、睡眠や体調の都合上実際のツイート間隔には比較的大きなばらつきがあった。しかし、50件のデータをもとに統計を取ったところ、確率変数「ツイート間隔の対数」が「一様分布に三次関数を適用したもの」として近似できること、その三次関数はツイート間隔を1時間の周辺に集中させるような性質を持つことに気づいた。以下、そのことについて解説と考察する。
    }
  >
  +section{データ}<
    +p{50件のデータは以下のとおりである。}
    +pn{0:59:42,
6:29:36,
0:54:56,
5:13:22,
1:13:42,
3:24:21,
1:04:07,
0:53:05,
0:59:23,
1:57:27,
1:14:50,
0:45:24,
0:59:06,
1:01:30,
8:54:04,
4:50:01,
0:58:24,
0:51:30,
0:49:54,
1:03:14,
1:28:24,
3:31:16,
3:06:19,
1:20:29,
1:36:08,
11:21:56,
0:41:41,
0:53:58,
2:25:50,
4:24:47,
1:48:54,
0:51:39,
8:33:02,
1:07:54,
1:06:54,
1:08:01,
1:05:13,
3:14:49,
3:01:49,
1:29:07,
3:06:35,
0:59:12,
1:46:46,
14:02:29,
0:52:34,
19:45:05,
1:07:07,
0:51:03,
1:12:10,
1:17:38,
}
    +pn{これを1時間で割って、ソートすると次のようになる。}
    +pn{
0.695,
0.757,
0.832,
0.851,
0.858,
0.861,
0.876,
0.885,
0.899,
0.916,
0.973,
0.985,
0.987,
0.990,
0.995,
1.025,
1.054,
1.069,
1.087,
1.115,
1.119,
1.132,
1.134,
1.203,
1.228,
1.247,
1.294,
1.341,
1.473,
1.485,
1.602,
1.779,
1.815,
1.958,
2.431,
3.030,
3.105,
3.110,
3.247,
3.406,
3.521,
4.413,
4.834,
5.223,
6.493,
8.551,
8.901,
11.366,
14.041,
19.751
    }
    +pn{ここで、このデータの自然対数が、 区間${ \[ 0 , 1 \] }上の連続一様分布に三次関数を適用して作られたものであるという仮定をおき、これら50個のデータの自然対数にそれぞれ${0.01, 0.03, ... 0.99}を対応させ、desmosでcurve-fittingをしたところ（ https://www.desmos.com/calculator/rcvlqujl9i ）、求める三次関数は${\app{g}{x} = -0.31942 + 1.97413x -4.66222x^2 + 5.94988x^3}であることが分かった。なお、${R^2=0.9939}であった。}
  >
  +section{三次関数そのものに関する考察}<
    +p{
      上記の三次関数は${x=0.26119}で変曲点を持つ。さらに、${x=0.26119}においてこの三次関数が取る値は${-0.01583}であり、非常に${0}に近い。
    }
    +pn{
      手動botを運用しているときに、「およそ1時間おきにツイートする」という意識のもとでツイートを行っていたことを考えると、この結果は自然である。単調増加する三次関数において最も変化率が小さくなるのは変曲点であり、これはつまりツイートの間隔の自然対数が変曲点の近くとなる確率が高いということを意味する。今回はツイート間隔を1時間で割ってから自然対数を取っているので、つまりツイート間隔が多くの場合1時間程度であるということを意味する。
    }
    +p{
      この考察をもとに、${\app{g}{x}}の形を${\app{g}{x} = a+bx+cx^2+dx^3}として4つのパラメータで書くのではなく${\app{g}{x} = A\paren{x-B}^3+C\paren{x-B}}としてパラメータを3つにしてみると（ https://www.desmos.com/calculator/nih5fwbofx ）、${A=6.17755, B=0.271117, C=0.764937}と表わせ、${R^2=0.9939}のままであった。
    }
  >
  +section{得られた確率変数に関する考察}<
    +p{
      この考察で得られた確率変数は、区間${ \[ 0 , 1 \] }上の連続一様分布に従う確率変数${K}を導入することで${X = e^{A\paren{K-B}^3+C\paren{K-B}}}と書ける。ただし、${A>0, C>0, 0<B<1}とする。
      ${B}は${X}が1以下になる確率を表す。${C}は何を表すだろうか？
    }
    +pn{
      ${K}が${B}に非常に近いときは${X \approx e^{C\paren{K-B}}}と考えることができるため、${X=1}周りでの累積頻度関数は${y=e^{C\paren{x-B}}}の逆関数${y=B+\frac{\ln x}{C}}で近似でき、故に確率密度関数はその導関数である${\app{f}{x} = \frac{1}{Cx}}と近似できる。したがって、${x=1}付近での確率密度関数の値の逆数が${C}であると考えることができる。（後に示すが、実は${x=1}においての確率密度関数の値は厳密に${\frac{1}{C}}である。）
    }
    +subsection{台}<
      +p{
        ${A\paren{x-B}^3+C\paren{x-B}}も${e^x}も単調増加するので、台は当然${\[e^{-AB^3-BC},e^{A\paren{1-B}^3+C\paren{1-B}} \] }である。今回のパラメータの場合、${\[0.719,19.100\]}である。
      }
    >
    +subsection{中央値、第1四分位点、第3四分位点}<
      +p{
        中央値、第1四分位点、第3四分位点は、それぞれ自明に${e^{A\paren{\frac{1}{2}-B}^3+C\paren{\frac{1}{2}-B}}}、${e^{A\paren{\frac{1}{4}-B}^3+C\paren{\frac{1}{4}-B}}}、そして${e^{A\paren{\frac{3}{4}-B}^3+C\paren{\frac{3}{4}-B}}}である。今回の場合、それぞれ${1.283}、${0.984}、${2.843}である。
      }
    >
    +subsection{期待値}<
      +p{
        区間${ \[ -B , 1-B \] }上の連続一様分布に従う確率変数${L}を導入する。
      }
      +math(${
        \app{E}{X} = \app{E}{
          e^{AL^3+CL}
        } = \int_{-B}^{1-B} e^{Ax^3+Cx} dx
        
      });
      +p{
        厳密解など出るはずもないので、近似する。
      }
      +math(${
        = \int_{-B}^{1-B} e^{Cx}\sums{k=0}{\infty} \frac{\paren{Ax^3}^k}{k\!} dx
        = \sums{k=0}{\infty}\frac{A^k}{k\!} \int_{-B}^{1-B} e^{Cx} x^{3k} dx
      });
      +pn{
        ここで${
          \int e^{Cx} x^{3k} dx = \frac{e^{Cx}}{C^{3k+1}}\sums{l=0}{3k}\paren{-1}^{3k-l}\paren{Cx}^{l}\frac{\paren{3k}\!}{l\!}
        }…はい、使い物になりませんね。おしまい。
      }
      +p{
        もうこれシンプソンの公式でいいんじゃないかな、上で第1四分位点とか出したんだし。${\app{G}{x} = e^{A\paren{x-B}^3+C\paren{x-B}}} とおけば、
      }
      +math(${
        \app{E}{X} = \int_{0}^{1} \app{G}{x} dx \approx \frac{1}{4 \cdot 3}\paren{\app{G}{0} + 4\app{G}{0.25} + 2\app{G}{0.5} + 4\app{G}{0.75} + \app{G}{1}}
        
      });
      +pn{今回のパラメータだと誤差は10\%ぐらい。十分許容範囲では？}
    >
    
  >
  +section{累積分布関数と確率密度関数}<
    +subsection{累積分布関数と確率密度関数を求める}<
      +p{
        ${\app{G}{x} = e^{A\paren{x-B}^3+C\paren{x-B}} \paren{0 < x < 1}}の逆関数が累積分布関数である。つまり、${x = e^{A\paren{y-B}^3+C\paren{y-B}}}なので
        ${
          \ln x = A\paren{y-B}^3+C\paren{y-B}
        } である。ここで
        ${
          k=\sqrt{\frac{3A}{4C}}
        } とおくと
      }
      +math(
        ${ 
          \frac{3k}{C}\ln x = 4k^3 \paren{y-B}^3+3k\paren{y-B}
        }
      );
      +pn{
        ここで${k\paren{y-B}=\sinh \lambda}とおくと${ 
          \frac{3k}{C}\ln x = 4\sinh^3 \lambda+3\sinh \lambda = \sinh 3\lambda
        }であり、これは${\frac{\app{\arsinh}{\frac{3k}{C}\ln x}}{3} = \lambda}と書ける。
      }
      +pn{
        ゆえに、累積分布関数は${y=\app{F}{x} = B + \frac{1}{k}\sinh \frac{\app{\arsinh}{\frac{3k}{C}\ln x}}{3}}(${\app{G}{0}<x<\app{G}{1}})であり、確率密度関数はこれを微分して
      }
      +math(${
        \app{f}{x} = \frac{1}{k}\frac{d}{dx}\sinh \frac{\app{\arsinh}{\frac{3k}{C}\ln x}}{3} = \frac{1}{k}\cosh \frac{\app{\arsinh}{\frac{3k}{C}\ln x}}{3} \cdot \frac{d}{dx} \frac{\app{\arsinh}{\frac{3k}{C}\ln x}}{3}
      });
      +math(${
       = \frac{1}{3k}\cosh \frac{\app{\arsinh}{\frac{3k}{C}\ln x}}{3} \cdot \frac{1}{\sqrt{1+\paren{\frac{3k}{C}\ln x}^2}} \cdot \frac{d}{dx} \frac{3k}{C}\ln x
      });
      +math(${
       = \frac{1}{Cx}\cosh \frac{\app{\arsinh}{\frac{3k}{C}\ln x}}{3} \cdot \frac{1}{\sqrt{1+\paren{\frac{3k}{C}\ln x}^2}}
      });
      +p{
        つまり、${\app{f}{x}}の式の形の中に露わには${B}が現れない。（${B}は定義域を定める役割のみを果たす。）
      }
      +pn{
        以後、${D=\frac{3k}{C}}なるパラメータを導入することとする。今回の例では${D=9.652}である。
      }
    >
    +subsection{パラメータと関数の分離}<
      +p{
        以上の結果の一部を https://www.desmos.com/calculator/lmw7kroctq にまとめた。図示すると、
        ${\paren{
          \app{G}{0},
          \app{f}{\app{G}{0}}
        }}から始まり、
        ${\paren{
          1,
          \app{f}{1}
        }}付近にあるピークまで下凸で登った後、一気に下凸で下っていく様子が見て取れる。この挙動を捉えたより単純な関数で確率密度関数を近似することを考える。
      }
      +pn{
        まず、${\app{H}{u} = \cosh \frac{\app{\arsinh}{u}}{3} \cdot \frac{1}{\sqrt{1+u^2}}}とおくと${\app{f}{x} = \frac{\app{H}{D\ln x}}{Cx}}と書ける。${\app{h}{x} =\app{H}{D\ln x}}とおくと${h}は${x=1}で極値${1}を取る。なぜなら、まず${\app{h}{1} =\cosh \frac{\app{\arsinh}{0}}{3} \cdot \frac{1}{\sqrt{1+0^2}}=1}であり、なおかつ

      }
      +math(${
        \frac{dh}{dx} = \frac{dH}{du}\cdot\frac{d\paren{D\ln x}}{dx} 
      });
      +math(${
        = \frac{D}{x}\paren{
          \paren{\frac{d}{du}\cosh \frac{\arsinh u}{3}} \cdot \frac{1}{\sqrt{1+u^2}} + 
          \cosh \frac{\arsinh u}{3} \cdot \frac{d}{du}\paren{\frac{1}{\sqrt{1+u^2}}}
        }
      });
      +math(${
        = \frac{D}{x}\paren{
           \frac{\app{\sinh}{\frac{1}{3} \arsinh u}}{ 3 \sqrt{u^2 + 1}} \cdot \frac{1}{\sqrt{1+u^2}} + 
          \cosh \frac{\arsinh u}{3} \cdot  -\frac{u}{\paren{1 + u^2}^{3/2}}
        }
      });
      +pn{
      は${u=D\ln x=0}で${0}になるからである。
      }
      +pn{
        さて、以上より、${\app{H}{u} = \cosh \frac{\app{\arsinh}{u}}{3} \cdot \frac{1}{\sqrt{1+u^2}}}とおくと${\app{f}{x} = \frac{\app{H}{D\ln x}}{Cx}}と書け、なおかつ${\app{H}{u}}は${u=0}で極値${1}を取ることが分かる。
      }
      +p{
        ここまでの結果をまとめ、${B, C, D, H}のみで書き直すと、次のようになる。
      }
      +listing{
        * ${\app{G}{x} = e^{\frac{4C^3D^2}{27}\paren{x-B}^3+C\paren{x-B}}} の逆関数 ${\app{F}{x} = B + \frac{3}{CD}\sinh \frac{\app{\arsinh}{D\ln x}}{3}}(${\app{G}{0}<x<\app{G}{1}})が累積分布関数である
        * 確率密度関数はその導関数である${\app{f}{x} = \frac{\app{H}{D\ln x}}{Cx}}(${\app{G}{0}<x<\app{G}{1}})
        * ただし、${\app{H}{u} = \cosh \frac{\app{\arsinh}{u}}{3} \cdot \frac{1}{\sqrt{1+u^2}}}
        * 特に、${\app{f}{1} = \frac{\app{H}{0}}{C} = \frac{1}{C}}であるので、パラメータ${C}は${\frac{1}{\app{f}{1}}}という意味合いを持った量である。
        * ${D}は、大きくすると${\app{G}{1}}を大きくさせ、また${\app{f}{x}}を${x=1}の辺りに尖らせる役割を持つ。
        * ${B}は、変わらず${X}が1以下になる確率を表す。これは${\app{F}{1} = B + \frac{3}{CD}\sinh \frac{\app{\arsinh}{0}}{3}=B}という式でも表現される。
      }
      +pn{
        これにより、パラメータの方は説明をつけることができたので、あとはパラメータを含まないこの謎の関数${\app{H}{u}}の性質を調べればよい。
      }
    >
  >
  +section{${\app{H}{u}}の性質}<
    +subsection{近似}<
      +pn{
        ${\app{\arsinh}{u} = \app{\ln}{u+\sqrt{u^2+1}}}だが、これを${\app{\ln}{2\sqrt{u^2+1}}}で近似しても${\app{H}{u}}の概形は${u=0}の周り以外ではほとんど変わらず、${u=0}の周りでも2\%ほどしか変わらない。}
      +pn{さて、その近似をした${\app{H'}{u}}は${\app{J}{v}=v^{2/3}+v^{4/3}}を用いて${\app{H'}{u} = \app{J}{\frac{1}{2\sqrt{1+u^2}}}}と書ける。何故かと言うと、
      }
      +math(${\app{H'}{u} = \frac{\app{\exp}{\frac{\app{\ln}{2\sqrt{u^2+1}}}{3}}+\app{\exp}{-\frac{\app{\ln}{2\sqrt{u^2+1}}}{3}}}{2}  \cdot \frac{1}{\sqrt{1+u^2}}});
      +math(${
        = \paren{\app{\exp}{\frac{\app{\ln}{v^{-1}}}{3}}+\app{\exp}{-\frac{\app{\ln}{v^{-1}}}{3}}}  \cdot v
      });
      +math(${
        = \paren{\app{\exp}{\app{\ln}{v^{-1/3}}}+\app{\exp}{\app{\ln}{v^{1/3}}}}  \cdot v = \app{J}{v}
      });
      +pn{
        と書けるからである。つまり、
      }
      +math(${
        \app{H}{u} \approx \paren{\frac{1}{4\paren{1+u^2}}}^{1/3}+\paren{\frac{1}{4\paren{1+u^2}}}^{2/3}
      });

    >
    +subsection{別方面}<
      +pn{
        ${u \to \infty}のとき${\app{\arsinh}{u} \approx \ln 2u}なので、${\cosh \frac{\app{\arsinh}{u}}{3} \approx \frac{u^{1/3}}{2^{2/3}}} である。故に、${\app{H}{u} \approx \paren{2u}^{-2/3}}である。したがって、${\paren{\app{H}{u}}^3 \approx \paren{2u}^{-2}}であることが分かる。
      }
      +pn{
      ${0}付近の挙動を調べるため、${\paren{\app{H}{u}}^{-3}}をテイラー展開すると、
      }
      +math(${
        \paren{\app{H}{u}}^{-3} = 1+\frac{4u^2}{3}+\frac{16u^4}{81}-\frac{128u^6}{2187}+\frac{512u^8}{19683} - ...
      });
      +pn{
        とまあ、なんとも興味深い形が現れる。（念のために言っておくと、この次の項でこのパターンは崩れ、2のべき乗でも3のべき乗でもない数が現れる。）
      }
      +pn{${v=\frac{4u^2}{9}}とすると、}
      +math(${
        \paren{\app{H}{u}}^{-3} = 1+3v+v^2-\frac{2}{3}v^3+\frac{2}{3}v^4 - ...
      });
      +pn{
        一方、これが${u \to \infty}のときは${\paren{\app{H}{u}}^{-3} \approx 4u^2}…とはいえ、これはあくまで漸近的な評価にすぎず、かなり悪い近似である。
      }
      +pn{
        もっと精密に評価しよう。${k=\paren{\frac{1}{4\paren{1+u^2}}}^{1/3}}としてやると
        ${
          \app{H}{u} \approx k+k^2
        }なので
      }
      +math(${
          \paren{\app{H}{u}}^{-3} \approx \paren{k^3+3k^4+3k^5+k^6}^{-1}
        }
      );
      +pn{
        だが、実際に近似してみると${k^6}を抜いたほうがむしろ精度が上がる気がする。（最初に「${\app{\ln}{2\sqrt{u^2+1}}}で近似」で導入された誤差を打ち消すように働く）
      }
      +pn{
        具体的に係数をどうすべきか考えてみよう。
      }
      +math(${
        \paren{\app{H}{u}}^{-3} = 1+3v+v^2-\frac{2}{3}v^3+\frac{2}{3}v^4 - ...
      });
      +pn{
        と
      }
      +math(${
          k^{-3}\paren{\app{H}{u}}^{3} \approx 1+3k+3k^2+k^3
        }
      );
      +pn{
        一方${k=\paren{\frac{1}{4+9v}}^{1/3}}より${k^{-3}=4+9v}である。前述のテイラー展開の逆数を取ると
      }
      +math(${
        \paren{\app{H}{u}}^{3} = 1-3v+8v^2-\frac{61}{3}v^3+\frac{151}{3}v^4 - ...
      });
      +pn{
        なので
      }
      +math(${
          k^{-3}\paren{\app{H}{u}}^{3} = 4 - 3v + 5v^2 - \frac{28}{3}v^3+\frac{55}{3}v^4 - ...
        }
      );
      +pn{これが}
      +math(${
          k^{-3}\paren{\app{H}{u}}^{3} \approx a+bk+ck^2+dk^3
        }
      );
      +pn{となるべく一致してほしい。まずは${d=0}で考える。${v=0}としたときの値を比較すると}
      +math(${
        4 = a+b\paren{\frac{1}{4}}^{1/3}+c\paren{\frac{1}{4}}^{2/3}
      });
      +pn{導関数の${v=0}のときの値を比較すると}
      +math(${
        -3 = b \frac{d\paren{4+9v}^{-1/3}}{dv} + c \frac{d\paren{4+9v}^{-2/3}}{dv} = -\frac{b}{3}\cdot 9\cdot \paren{4+9v}^{-4/3}-\frac{2c}{3}\cdot 9\cdot \paren{4+9v}^{-5/3}
      });
      +pn{一旦${c=0}とすると${
        1 = b \cdot \paren{4}^{-4/3}
        }と${4 = a+b\paren{\frac{1}{4}}^{1/3}}より${b=4^{4/3}, a = 4 -4^{4/3}\paren{\frac{1}{4}}^{1/3}= 0}、ふーむ。
      }
      +pn{グラフ書いてみたんですが、えーなんかめっちゃ合いますね、完璧では。}
      +pn{ 
        ${
          \app{H}{u} \approx 4^{4/9}k^{4/3} = \paren{\frac{4}{4+9v}}^{4/9} = \paren{\frac{1}{1+u^2}}^{4/9}
        }
      }
      +pn{
        1から遠いところではズレないわけではないけどまあ仕方がない。
      }
      +pn{
        2階微分まで使うとか、漸近評価も使うとかすればいいんだろうけど。
      }
    >
  >
 +section{近似した${\app{H}{u}}で考察}<
  +subsection{累積分布関数}<
  +p{ 
    確率密度関数が${\app{f}{x} = \frac{\app{H}{D\ln x}}{Cx}}、${\app{H}{u} = \paren{\frac{1}{1+u^2}}^{4/9}} である場合、
  }
  +math(  ${
      \app{F}{x'} =\int_{a}^{x'} \frac{\app{H}{D\ln x}}{Cx}dx = \int_{D\ln a}^{D \ln x'} \frac{\app{H}{q}}{CD}dq = \int_{D\ln a}^{D \ln x'} \frac{\paren{1+q^2}^{-4/9}}{CD}dq
    });
  +pn{
    …初等関数で積分できないなこれ。一体なんのための近似だったのやら。
  }
  >
  +subsection{期待値}<
    +p{
      期待値なら${x}が入るから積分できるのでは？
    }
    +pn{
      と一瞬期待してしまったが、${\ln x}が入るので普通に無理である。
    }
  >
 >
 +section{そもそも三次関数を使ったのが良くなかったのでは？}<
 >
>