Qt案外使いやすいかも

ずっと敬遠していたGUIツールキットQt。

敬遠してたのは開発環境やら何やらごてごてしてるから。

もっと気軽にlibファイルだけリンクして使えるシンプルなものが好きなんだが。

でもちょっと使ってみると案外使いやすいかも。

QtCreatorも結構いい感じかも。

 

ちょっとWindows環境でGUIプログラム作る必要が生じていじっている。

QtCreator、ちょっと慣れるのに時間がかかったけど案外Visual Studioより使いやすいかも。

Linuxでも使ってみようかな。

ジャジャガッチ | C/C++ | 22:14 | comments(0) | trackbacks(0) |

小林泰三「セピア色の凄惨」読了

小林泰三「セピア色の凄惨」読了。

やっぱり小林泰三、安定して面白い。

ちょっと(かなり?)病んだり狂ったりした人たちの短編集。

気軽に読める長さだけど内容は暗くて面白い。

ジャジャガッチ | その他 | 22:09 | comments(0) | trackbacks(0) |

フィッティングで方程式を解いてみる(3)デノイズ

前回、前々回、フィッティングで方程式を解いてみた。

例が簡単過ぎてあまり面白みを感じない方もいそうなので、もう少し凝ったことをしてみたい。

 

今回は、記事タイトルに偽りありであるが関数の最小化をしてみたい。

フィッティングが二乗和で表される関数の最小化をしていることを利用して方程式を解いていたのだが、今回は方程式ではなく、関数の最小化問題を解いてみたいと思う。

 

今回行うのは簡単なデノイズである。

いま、あるデータ{Xi}がある。本来、データは滑らかでなければならないとする。この仮定の下でノイズを含まないデータの推定を行ってみよう。

デノイズ後の値を{ai}とする。滑らかさの指標として、

¥[
¥sum_i(a_{i+1}-a_i)^2
¥]

を用いる。これが小さければ滑らかである。最も滑らかなのは全ての点で同一の値を持つ場合である。これではデノイズ出来たとはいい難い。{ai}は{Xi}にそれなりに近い必要がある。そこで

¥[
¥sum_i(a_i-X_i)^2 + ¥lambda¥sum_i(a_{i+1}-a_i)^2
¥]

を最小化する{ai}を求めることにする。λは定数で、これを大きくすると、滑らかであるが元のデータとの乖離が大きくなる。

この関数は二乗和となっているので、フィッティング機能を利用して解くことが出来る。

 

フィッティング対象のデータは次のようなものを用意する。データの次元をNとする。

 

0 X[0]

1 X[1]

・・・・

N-1 X[N-1]

N 0

・・・

2N-2 0

 

最初のN行は

¥[
¥sum_i(a_i-X_i)^2 
¥]

の部分である。後半は

¥[
¥lambda¥sum_i(a_{i+1}-a_i)^2
¥]

の部分である。

 

フィッティング関数は

¥[
¥sum_{i=0}^{N-1}¥delta_{xi}a_i + ¥sqrt{¥lambda}¥sum_{i=N}^{2N-2}¥delta_{xi}(a_{i-N+1}-a_{i-N})
¥]

と定義する。これで元の問題をフィッティング問題として定義できた。

 

それではやってみよう。今回は問題が少し複雑で流石にGnuplotに手打ちでは辛いので、次のコードでテキストを生成してGnuplotで読み込んだ。

 

#define _CRT_SECURE_NO_WARNINGS
#include <stdio.h>
#include <vector>
#include <random>
#include <cmath>
double M_PI = 3.14159265358979;
using namespace std;

int main()
{

    FILE *fp = fopen("data.txt", "w");
    int N = 64;
    vector<double> val(N);
    std::mt19937 mt(0);
    std::normal_distribution<double> gauss(0, 0.1);
    for (int n = 0; n < N; n++)
    {
        val[n] = sin(2 * M_PI*n / (double)N) +gauss(mt);
        fprintf(fp, "%d %e¥n",n,val[n]);
    }
    for (int n = N; n < 2*N-1; n++)
    {
        fprintf(fp, "%d 0¥n",n);
    }
    fclose(fp);

    fp = fopen("macro.txt", "w");
    fprintf(fp,"d(i,j) = (i==j) ? 1 : 0¥n");
    fprintf(fp,"fit ");
    for (int n = 0; n < N; n++)
    {
        fprintf(fp,"d(x,%d)*a%d + ",n,n);
    }
    double lambda = 10;
    for (int n = N; n < 2*N-1; n++)
    {
        fprintf(fp,"%f*d(x,%d)*(a%d-a%d)", sqrt(lambda),n, n+1-N,n-N);
        if (n != 2 * N - 2) fprintf(fp," + ");
        else fprintf(fp," ");
    }
    fprintf(fp,"¥"data.txt¥" via ");
    for (int n = 0; n < N; n++)
    {

        fprintf(fp,"a%d ",n);
        if (n != N - 1) fprintf(fp,", ");

    }
    fclose(fp);

    return 0;
}

 

このコードは、オリジナルデータとしてsin関数を用意し、σ=0.1のガウスノイズを足し合わせてデータを生成している。

その後、Gnuplot用のコマンドを生成している。

具体的にはN=64, λ=1では次のコマンドを生成している。

d(i,j) = (i==j) ? 1 : 0
fit d(x,0)*a0 + d(x,1)*a1 + d(x,2)*a2 + d(x,3)*a3 + d(x,4)*a4 + d(x,5)*a5 + d(x,6)*a6 + d(x,7)*a7 + d(x,8)*a8 + d(x,9)*a9 + d(x,10)*a10 + d(x,11)*a11 + d(x,12)*a12 + d(x,13)*a13 + d(x,14)*a14 + d(x,15)*a15 + d(x,16)*a16 + d(x,17)*a17 + d(x,18)*a18 + d(x,19)*a19 + d(x,20)*a20 + d(x,21)*a21 + d(x,22)*a22 + d(x,23)*a23 + d(x,24)*a24 + d(x,25)*a25 + d(x,26)*a26 + d(x,27)*a27 + d(x,28)*a28 + d(x,29)*a29 + d(x,30)*a30 + d(x,31)*a31 + d(x,32)*a32 + d(x,33)*a33 + d(x,34)*a34 + d(x,35)*a35 + d(x,36)*a36 + d(x,37)*a37 + d(x,38)*a38 + d(x,39)*a39 + d(x,40)*a40 + d(x,41)*a41 + d(x,42)*a42 + d(x,43)*a43 + d(x,44)*a44 + d(x,45)*a45 + d(x,46)*a46 + d(x,47)*a47 + d(x,48)*a48 + d(x,49)*a49 + d(x,50)*a50 + d(x,51)*a51 + d(x,52)*a52 + d(x,53)*a53 + d(x,54)*a54 + d(x,55)*a55 + d(x,56)*a56 + d(x,57)*a57 + d(x,58)*a58 + d(x,59)*a59 + d(x,60)*a60 + d(x,61)*a61 + d(x,62)*a62 + d(x,63)*a63 + 1.000000*d(x,64)*(a1-a0) + 1.000000*d(x,65)*(a2-a1) + 1.000000*d(x,66)*(a3-a2) + 1.000000*d(x,67)*(a4-a3) + 1.000000*d(x,68)*(a5-a4) + 1.000000*d(x,69)*(a6-a5) + 1.000000*d(x,70)*(a7-a6) + 1.000000*d(x,71)*(a8-a7) + 1.000000*d(x,72)*(a9-a8) + 1.000000*d(x,73)*(a10-a9) + 1.000000*d(x,74)*(a11-a10) + 1.000000*d(x,75)*(a12-a11) + 1.000000*d(x,76)*(a13-a12) + 1.000000*d(x,77)*(a14-a13) + 1.000000*d(x,78)*(a15-a14) + 1.000000*d(x,79)*(a16-a15) + 1.000000*d(x,80)*(a17-a16) + 1.000000*d(x,81)*(a18-a17) + 1.000000*d(x,82)*(a19-a18) + 1.000000*d(x,83)*(a20-a19) + 1.000000*d(x,84)*(a21-a20) + 1.000000*d(x,85)*(a22-a21) + 1.000000*d(x,86)*(a23-a22) + 1.000000*d(x,87)*(a24-a23) + 1.000000*d(x,88)*(a25-a24) + 1.000000*d(x,89)*(a26-a25) + 1.000000*d(x,90)*(a27-a26) + 1.000000*d(x,91)*(a28-a27) + 1.000000*d(x,92)*(a29-a28) + 1.000000*d(x,93)*(a30-a29) + 1.000000*d(x,94)*(a31-a30) + 1.000000*d(x,95)*(a32-a31) + 1.000000*d(x,96)*(a33-a32) + 1.000000*d(x,97)*(a34-a33) + 1.000000*d(x,98)*(a35-a34) + 1.000000*d(x,99)*(a36-a35) + 1.000000*d(x,100)*(a37-a36) + 1.000000*d(x,101)*(a38-a37) + 1.000000*d(x,102)*(a39-a38) + 1.000000*d(x,103)*(a40-a39) + 1.000000*d(x,104)*(a41-a40) + 1.000000*d(x,105)*(a42-a41) + 1.000000*d(x,106)*(a43-a42) + 1.000000*d(x,107)*(a44-a43) + 1.000000*d(x,108)*(a45-a44) + 1.000000*d(x,109)*(a46-a45) + 1.000000*d(x,110)*(a47-a46) + 1.000000*d(x,111)*(a48-a47) + 1.000000*d(x,112)*(a49-a48) + 1.000000*d(x,113)*(a50-a49) + 1.000000*d(x,114)*(a51-a50) + 1.000000*d(x,115)*(a52-a51) + 1.000000*d(x,116)*(a53-a52) + 1.000000*d(x,117)*(a54-a53) + 1.000000*d(x,118)*(a55-a54) + 1.000000*d(x,119)*(a56-a55) + 1.000000*d(x,120)*(a57-a56) + 1.000000*d(x,121)*(a58-a57) + 1.000000*d(x,122)*(a59-a58) + 1.000000*d(x,123)*(a60-a59) + 1.000000*d(x,124)*(a61-a60) + 1.000000*d(x,125)*(a62-a61) + 1.000000*d(x,126)*(a63-a62) "data.txt" via a0 , a1 , a2 , a3 , a4 , a5 , a6 , a7 , a8 , a9 , a10 , a11 , a12 , a13 , a14 , a15 , a16 , a17 , a18 , a19 , a20 , a21 , a22 , a23 , a24 , a25 , a26 , a27 , a28 , a29 , a30 , a31 , a32 , a33 , a34 , a35 , a36 , a37 , a38 , a39 , a40 , a41 , a42 , a43 , a44 , a45 , a46 , a47 , a48 , a49 , a50 , a51 , a52 , a53 , a54 , a55 , a56 , a57 , a58 , a59 , a60 , a61 , a62 , a63 

 

このコマンドを実行すればGnuplotがフィッティングを実行する。

さて、結果を見てみよう。

まず、最初は理想データとノイズ入りデータである。

ガタガタしている。

λ=1でやってみよう。

少し滑らかになったのがわかるだろうか。

さらにλ=10ならばより滑らかになる。

目的関数の値もきちんと減少していることを確認した。

 

どうだろう、一応それなりに実践的なことも出来ることを示せたわけだが。面白くない?

役に立つかというと役には立たないわけだけど。

コード書かなくても最適化問題が解けるわけだが、単純な例なら今はWolfram alpha等便利なものがいくらでもあるし、複雑なものになると結局コードでコマンド生成しないと辛いものが・・・。

でも面白いからいいんだ。

ジャジャガッチ | 数学 | 22:42 | comments(0) | trackbacks(0) |

フィッティングで方程式を解いてみる(2)連立方程式

前回、f(x)=cをフィッティングで解いてみたわけだが、今回は連立方程式でやってみる。

基本は前回と同じである。

解きたい方程式を

¥[
¥sum_j A_{ij}x_j = b_i
¥]

とする。

次のようなデータを用意する。

0 b0

1 b1

2 b2

・・・

 

フィッティング関数は

¥[
f(x;¥{a¥}) = ¥sum_j A_{ij}a_j¥delta_{xi}
¥]

とする。aが推定するパラメータ、δはクロネッカーのデルタである。

フィッティングすれば

¥[
 ¥sum_i(¥sum_j A_{ij}a_j-b_i)^2
¥]

を最小とする解が求まる。これにより元の方程式の解が求まる。

具体的にやってみよう。

x+y=3

2x-y=0

を解いてみる。明らかに解はx=1, y=2である。

”data.txt”として

0 3
1 0

を用意する。

クロネッカーのデルタを定義する。

d(i,j) = (i==j) ? 1 : 0 

後は上で述べたようにフィッティング関数を定義してフィッティングすればよい。

 fit d(x,0)*a+d(x,0)*b+2*d(x,1)*a-d(x,1)*b "data.txt" via a,b

はい、求まった。

ジャジャガッチ | 数学 | 23:17 | comments(0) | trackbacks(0) |

フィッティングで方程式を解いてみる

実に半年ぶりのブログ更新。

Fallout4やったり、ウィッチャー3やったりしてた。

あんまり家でガッツリ趣味やる時間がなくて放置していた。

会社ではいろいろ面白いことやってるんだけどそれをここに書くわけにはいかないし。

趣味やるくらいなら英語やらねばというのもある。ここしばらく英語は勉強してなくて、3月に入ってから勉強を再開した。

 

本題。

今日は以前ふと思いついたことを実際にやってみることにした。

それは・・・、フィッティングを使って方程式を解くこと!

さて、どういうことをするかは後で説明することにして思い出話から。

高校の数学では、基本的に解析的に解けない方程式は出てこない。

しかし、大学に入ると解析的に解けない方程式のオンパレード。

ここで面食らうわけだが、大学1年生のときに先生にいいことを教わった。

グラフソフトでグラフ書けば数値解簡単に求められるよ、と。

例えば

x^2 -2x = -1

を考えよう。

この数値解はグラフソフトでy=x^2-2xとy=-1をプロットして交点の座標を読み取ればいい。

この例は解析的に解けるが、解析的に解けない方程式も同じだ。

 

さて、ここからが本題。

同じグラフソフトを使って解を求めるならフィッティング機能使えないかな?

グラフソフトにはフィッティング機能がついていることがある。

例えば僕が愛用しているGnuplotは

fit a*x "data.txt" via a 

と入力すれば直線フィットが出来る。

与えたデータに最も近くなるようにパラメータを推定してくれるわけだ。

これを利用して方程式を解けないだろうか。

 

フィッティングは次のようなことをしている。

対象のデータを{xi,yi}としよう。これをパラメータ{a}の関数fでフィッティングするとは

¥[
¥sum_i (y_i-f(x_i;¥{a¥}))^2
¥]

を最小化する{a}を求めることに他ならない。

 

それではこれを利用して方程式を解いてみよう。

ここでは

g(a) = c

をaについて解くことにする。cは定数である。

Gnuplotではデータは2点以上必要である。

次のようにデータを用意する。

0 0

1 c

フィッティング関数は

f(x;a) = g(a)*x

とする。

するとフィッティングすることにより求まるパラメータaは

¥[
(g(a)¥times 0-0)^2 + (g(a)-c)^2
¥]

を最小化するものとなる。これはg(a)=cである。

つまり数値解が求まる。

 

さて、x^2-2*x = -1をこの方法で解いてみよう。

"data.txt"として

0 0
1 -1

を用意する。

fit (a**2-2*a)*x "data.txt" via a

を実行すれば

はい、ばっちしa=1が得られました!

いや、だから何?なんだけど(笑)

役に立たなくても面白いものは面白いじゃない?

ジャジャガッチ | 数学 | 22:39 | comments(0) | trackbacks(0) |
1/1PAGES | |

03
--
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
--
>>
<<
--
PR
RECOMMEND
RECENT COMMENT
MOBILE
qrcode
OTHERS
Since 2013/09/17
LATEST ENTRY
CATEGORY
ARCHIVE
LINKS
PROFILE
SEARCH