memo

Transcription

1 数理情報工学演習第一 C ( 第 12 回 ) 2017/07/11 DEPARTMENT OF MATHEMATICAL INFORMATICS 1

2 今日の内容 : ファイルの入出力 コマンドライン引数 2 分探索 最長単調増加列 2

3 ファイル操作の手続き : ファイル操作 ファイルからのデータ読み込み ファイルへのデータ書き出し 基本的な手順 読みこむ / 書き出すファイルを開く (fopen) 読み込み / 書き出しの操作 (fscanf/fprintf) ファイルを閉じる (fclose) 3

4 ファイルを開く / 閉じる : fopen/fclose fopen( ファイル名, モード ) モード :r( 読み込み ), w( 書き込み ), a( 追記 ) FILE *fp; // ファイルポインタの宣言 fp = fopen( input.txt, r ) ; // 読み込みモードで input.txt を開く 各種ファイル処理 fclose(fp); // 閉じる 4

5 ファイルへの書き込み 読み込み : ファイルポインタ以外はprintf/scanfと同様 書き出し :fprintf( ファイルポインタ, ) fprintf(fp, %d n, i) 読み込み :fscanf( ファイルポインタ, ) 整数は %d, 実数は %lf int x; double y; fscanf(fp, %d %lf, &x, &y); ファイルの最後に到達した場合 EOF を返す if ( fscanf( ) == EOF ){.} 5

6 #include <stdio.h> #include <stdlib.h> void read_file(char *filename) { double x; FILE *fp; fp = fopen(filename, r ); if (fp == NULL) { printf( cannot open %s n, filename); exit(1); } while (!feof(fp)) { if (scanf("%lf", &x) < 1) { break; } fprintf(stdout, %lf n, x); } fclose(fp); } int main(argc, char *argv[]) { read_file( test.txt ); return 0; } 6

7 コマンドライン引数の利用./a.out input.txt と引数をつけて実行したとき argc = 4; argv[0] =./a.out ; プログラム名 argv[1] = input.txt ; argv[2] = 10 ; argv[3] = 0.1 ; int main(int argc, char *argv[]) { FILE *f; int n; double x; if (argc < 4) exit(1); // 引数が足りない f = fopen(argv[1], r ); n = atoi(argv[2]); // 文字列を整数に x = atof(argv[3]); // 文字列を実数に } 7

8 二分探索 アルゴリズムとデータ構造で重要な概念 全順序集合の探索を高速化する 集合 S の要素を L, E, G に分ける L = {x x S, x < p} (p より小さい要素 ) E = {x x S, x = p} (p と等しい要素 ) G = {x x S, x > p} (p より大きい要素 ) k を探索するとき p = k ならば探索終了 ( 見つかった ) p < k ならば G を二分探索 p > k ならば L を二分探索 8

9 既ソート配列での二分探索 E は配列の中央の要素 (p = S[n/2]) L は中央より左側の要素 (S[0..n/2 1]) G は中央より右側の要素 (S[n/2+1..n 1]) 集合が配列のl 番目からh 番目で表されているとき m = (l+h) / 2 とする S[m]=k ならば探索終了 (k が存在した ) S[m]<k ならば k はL, E には存在しない G を探索 S[m]>k ならば k はG, E には存在しない L を探索 l m h < p p > p 9 L E G

10 int bsearch(int n, int *S, int k) { int high, low, mid; low = 0; // 探索の範囲は最初は配列全体 [0..n-1] high = n-1; while (low <= high) { mid = (low + high) / 2; if (S[mid] == k) { // 見つかった return mid; } else if (S[mid] < k) { low = mid + 1; // 探索範囲を後半だけにする } else { high = mid - 1; // 探索範囲を前半だけにする } } return -1; // 見つからなかった } 計算時間 O(log n) 10

11 このアルゴリズムでは,S[i] = k となる i が複数ある場合にどれを返すか分からない. S[0,n-1] で,S[i] >= k となる最小の i を求めることを考える そのような i が無いときは n を返すとする S[n] = とみなしていることになる While ループ実行中に, 常に次の条件を満たすようにする i の探索範囲は [low, high] 探索範囲内に答えが無い場合, 答えは範囲のすぐ右にある 最初の探索範囲が [0, n-1] の時, この条件を満たす 探索範囲を狭めるときにこの条件が成り立つことを確認 探索範囲が空集合になった時に何を返せばいいか考える 11

12 次のようなプログラムで動作確認をする int main(int argc, char *argv[]) { int i, pos, k, n; int *X; FILE *fin; fin = fopen(argv[1], "rb"); fscanf(fin, " %d", &n); X = (int *)malloc(n * sizeof(int)); for (i=0; i<n; i++) { fscanf(fin, " %d", &X[i]); } k = atoi(argv[2]); pos = binary_search(n, X, k); for (i=0; i<n; i++) { printf(" %d", X[i]); } printf(" n"); printf("key = %d pos = %d n", k, pos); } test.txt 実行結果 /a.out test.txt key = -1 pos = 0./a.out test.txt key = 2 pos = 3./a.out test.txt key = 6 pos = 6./a.out test.txt key = 7 pos = 6./a.out test.txt key = 8 pos = 9./a.out test.txt key = 10 pos = 11 12

13 最長単調増加列 入力 : 整数の配列 A[1,n] 出力 : A[i1] < A[i2] < < A[ik] (1 i1 < i2 < < ik n) となるような最大の k つまり,A の部分列で狭義単調増加になっているもので最長のものの長さ 単純にやると, 答えの候補は 2 n 個 動的計画法を用いる 13

14 動的計画法による解法 1 長さ n+1 の配列 L と長さ n の配列 P を用いる L[i] = A[1] から A[i] の部分列で, 最後が A[i] である解の長さ P[i] = その解で A[i] の1つ前に現れる値 A[p] の添え字 p L[0] = 0; A[0] = とする L[i] の最大値が答え ( 長さ 4, 部分列 1, 2, 5, 7) A L P

15 1,,i-1 まで計算してある状態で,P[i] と L[i] を求めるには 0 j i-1 に対し,A[j] が最後の要素である最長単調増加列の後に A[i] を追加する. つまり A[j] < A[i] となる j で L[j] が最大のものの後に A[i] を追加する. L[i] = L[j] + 1 P[i] = j 計算時間は O(n 2 ) 15

16 A L P i = 6 のとき,A[6] = 5 未満の A[j] (0 j < i) は A[0], A[1], A[2], A[4] その中で L が最大の A[4] の後に A[6] をつけるつまり P[6] = 4 L[6] = L[4]

17 動的計画法による解法 2 長さ n+1 の配列 M と長さ n の配列 P を用いる これまでに見つかった最長単調増加列の長さを len とする A[i] までの処理を行った段階で M[j] = A[0] から A[i] の長さが j の部分列で最後の値 A[k] が最小である k P[k] = その解で A[k] の1つ前に現れる値 A[p] の添え字 p M[0] = -1 とする i=0 A M P

18 M[j] = 長さが j の部分列で最後の値 A[k] が最小である k P[k] = その解で A[k] の 1 つ前に現れる値 A[p] の添え字 p i=1 A M P i=2 A M P i=3 A M P

19 M[j] = 長さが j の部分列で最後の値 A[k] が最小である k P[k] = その解で A[k] の 1 つ前に現れる値 A[p] の添え字 p i=4 A M P i=5 A M P i=6 A M P

20 M[j] = 長さが j の部分列で最後の値 A[k] が最小である k P[k] = その解で A[k] の 1 つ前に現れる値 A[p] の添え字 p i=7 A M P 長さが 4: A[6] = 7, A[5] = 5, A[3] = 2, A[1] = 1 長さが 3: A[7] = 3, A[3] = 2, A[1] = 1 長さが 2: A[3] = 2, A[1] = 1 長さが 1: A[1] = 1 A[M[j]] (1 j len) は単調増加 A[M[1]] = 1, A[M[2]] = 2, A[M[3]] = 3, A[M[4]] = 7 なぜならば, もし長さ j 2 の A[M[j]] で終わる単調増加列があるなら, A[P[M[j]]] で終わる長さ j-1 の単調増加列があり, A[M[j]] > A[P[M[j]]] A[M[j-1]] 20

21 つまり,A[i] を追加する際, A[M[j-1]] < A[i] A[M[j]] となる j を求めるそのような 1 j len がなければ,j = len+1 とする M[j] = i とする j > len のときは len を更新 (len = j) 計算時間は O(n log n) 21

22 i=2 A M P A[M[j]] 1 6 (len = 2) A[3] = 2 A[M[2]] = 6 なので, j = 2 i= A M P A[M[j]] 1 2 (len = 2) 22

23 i=3 A M P A[M[j]] 1 2 (len = 2) A[4] = 8 > A[M[2]] = 2 なので, j = 3 i= A M P A[M[j]] (len = 3) 23

24 課題 最長単調増加列の長さを求める 入力 : 整数列 ( 標準入力より ) 1 行目 : 整数の個数 n 2 行目 : n 個の整数 出力 : 最長単調増加列の長さ増加列自体は出力しない まず解法 1 で解き, 時間があれば解法 2 で再提出 7/11( 火 ) 17:55 までに提出 24

25 文字列の編集距離 文字列 X と Y の編集距離 (edit distance) とは, X に以下の編集操作を繰り返して Y に変換する際の最小の操作回数である. 挿入 : x i と x i+1 の間に文字 c を入れる 削除 : x i を削除 置換 : x i を文字 c で置き換える X = ACGTT CGTT CGCTT Y = CGCAT A を削除 C を挿入 T を A に置換 編集距離 = 3 25

26 X の文字を削除する代わりに,Y に gap ( 空白 ) を入れる X に文字を挿入するときは必ず Y の文字を入れることになる X に gap を入れる X と Y に gap を入れた X と Y で, ミスマッチの数を最小にする問題と等価 X = ACGTT Y = CGCAT X = ACG TT Y = CGCAT ミスマッチ = 3 26

27 X i = <x 1, x 2,, x i > と Y j = <y 1, y 2,, y j > の編集距離を c[i,j] とする x i と y j をマッチさせる場合 x i = y j ならば c[i,j] = c[i 1,j 1] x i y j ならば c[i,j] = c[i 1,j 1]+1 x i の次に gap を入れ y j とマッチさせる場合 c[i,j] = c[i,j 1]+1 y j の次に gap を入れ x i とマッチさせる場合 c[i,j] = c[i 1,j]+1 c[i,j] はこれら3つの中の最小値 O(mn) 時間 (m = X, n= Y ) 27

28 28 は Y に gap を入れることを表す は X に gap を入れることを表す 斜めは一致または置換を表す x i A C G T T i y j C G C A T j X = ACGTT Y = CGCAT X = ACG TT Y = CGCAT

29 実際にはミスマッチの数ではなく, 以下のようなコストを考える x i と y j をマッチさせる場合 x i = y j ならば c[i,j] = c[i 1,j 1] x i y j ならば c[i,j] = c[i 1,j 1] +1 x i の次に gap を入れ y j とマッチさせる場合 c[i,j] = c[i,j 1]+2 y j の次に gap を入れ x i とマッチさせる場合 c[i,j] = c[i 1,j]+2 注 : このアルゴリズムは任意のアルファベットで動作する 29

30 実は, 次のグラフで (0,0) から (n, m) への最短路長を求めればよい j i y j C G C A T x i A C G A A T 2 4 T 5 T 30

31 ヒント 次のような関数を作成する graph *graph_new(int n); 節点数 n で, 枝の無いグラフを作成する 文字列長を n, m とすると, 節点数は (n+1)(m+1) なので注意 void graph_insert_edge(graph *G, int u, int v, double c); グラフに, 有向枝 (u, v) を追加. 枝の長さは c これらの関数を使い, 格子状のグラフを作る 31

32 ヒント 文字列は次のようにして読み込む int read_file(char *fname, char **buf) // ファイル名 fname, 読み込む配列 buf ( この関数内でメモリ確保 ) { int n; char *p; FILE *fin; fin = fopen(fname, "r"); // ファイルを読み込みモードでオープン if (fin == NULL) { // オープンできなければエラー終了 printf("cannot open %s n", fname); exit(1); } fseek(fin, 0, SEEK_END); n = ftell(fin); // ファイルの長さを取得 fseek(fin, 0, SEEK_SET); p = malloc(n*sizeof(char)); // ファイルを格納するメモリを確保 if (p == NULL) { printf("not enough memory n"); exit(1); } fread(p, sizeof(char), n, fin); // メモリ p に n 文字読み込む fclose(fin); // ファイルを閉じる *buf = p; // 文字列を格納するポインタを関数の呼び出し側に返す return n; // 文字列長を返す } 32

33 int main(int argc, char *argv[]) { int n, m; int i, c; FILE *fin; char *X, *Y; if (argc < 3) { // コマンドライン引数が足りなければ終了 printf("%s file1 file2 n", argv[0]); // argv[0] はプログラム名が入っている return 1; } n = read_file(argv[1], &X); // 1 つ目のファイルを X に読み込む m = read_file(argv[2], &Y); // 2 つ目のファイルを Y に読み込む c = sm(n, m, X, Y); // X と Y の編集距離を計算 free(x); free(y); printf("%d n", c); // 答えを出力 return 0; } 33

34 最終レポート 2つの文字列間の編集距離を求め, 実行時間を計測する 動的計画法と最短路に基づく手法の速度比較をする 異なる長さの文字列を複数用い, 実行時間の増加を調べる 用いる文字列は自分で用意する DNAでもそうでなくても良い 締切 :7/24( 月 ) 宛先 :miprogramming2017+final@gmail.com 34

35 レポートについての注意 メールの題名 : 学生証番号 8 桁ハイフンなし ( 半角スペース ) 氏名 ( 半角スペース ) 最終レポート提出 例 : 定兼邦彦最終レポート提出 添付ファイル : レポート ( 学生証番号 8 桁.pdf 5 枚以内 ) ソースコードを圧縮したもの ( 学生証番号 8 桁.zip) 例 : pdf 例 : zip レポートについて タイトル 所属 学生証番号 氏名 日付を一ページ目に書くこと ( 独立させたページにはしなくてよい 1. 課題内容 2. 手法 ( データ構造 アルゴリズム ) の説明 3. 結果 4. 考察の章立てで書くこと A4 5 枚以下に収めること 以下は本レポートに限らず留意のこと グラフにはタイトル 縦軸 横軸 ( 単位も ) 凡例を忘れずに 表のタイトルを忘れずに グラフ 表のタイトルは本文を読まなくても何が書いてあるかわかるように 35