Context Navigation

main.cpp @ 1343

Revision 1343, 13.2 kB (checked in by sindj, 13 years ago)
Importance sampling. JS

Line
1
2	/*!
3	\file
4	\brief Robust
5	\author Vasek Smidl
6
7	*/
8
9	#include "estim/arx.h"
10	#include "robustlib.h"
11	#include <vector>
12	#include <iostream>
13	#include <fstream>
14	#include <itpp/itsignal.h>
15
16	using namespace itpp;
17	using namespace bdm;
18
19	const int emlig_size = 2;
20	const int utility_constant = 2;
21
22
23	int main ( int argc, char* argv[] ) {
24
25	itpp::Laplace_RNG LapRNG = Laplace_RNG();
26
27	/*
28	// EXPERIMENT: 100 AR model generated time series of length of 30 from y_t=0.95y_(t-1)+0.05y_(t-2)+0.2*e_t,
29	// where e_t is normally, student(4) and cauchy distributed are tested using robust AR model, to obtain the
30	// variance of location parameter estimators and compare it to the classical setup.
31	vector<vector<vector<string>>> string_lists;
32	string_lists.push_back(vector<vector<string>>());
33	string_lists.push_back(vector<vector<string>>());
34	string_lists.push_back(vector<vector<string>>());
35
36	char* file_strings[3] = {"c:\\ar_normal.txt", "c:\\ar_student.txt", "c:\\ar_cauchy.txt"};
37
38
39	for(int i = 0;i<3;i++)
40	{
41	ifstream myfile(file_strings[i]);
42	if (myfile.is_open())
43	{
44	while ( myfile.good() )
45	{
46	string line;
47	getline(myfile,line);
48
49	vector<string> parsed_line;
50	while(line.find(',') != string::npos)
51	{
52	int loc = line.find(',');
53	parsed_line.push_back(line.substr(0,loc));
54	line.erase(0,loc+1);
55	}
56
57	string_lists[i].push_back(parsed_line);
58	}
59	myfile.close();
60	}
61	}
62
63	for(int j = 0;j<string_lists.size();j++)
64	{
65
66	for(int i = 0;i<string_lists[j].size()-1;i++)
67	{
68	vector<vec> conditions;
69	//emlig* emliga = new emlig(2);
70	RARX* my_rarx = new RARX(2,30);
71
72	for(int k = 1;k<string_lists[j][i].size();k++)
73	{
74	vec condition;
75	//condition.ins(0,1);
76	condition.ins(0,string_lists[j][i][k]);
77	conditions.push_back(condition);
78
79	//cout << "orig:" << condition << endl;
80
81	if(conditions.size()>1)
82	{
83	conditions[k-2].ins(0,string_lists[j][i][k]);
84
85	}
86
87	if(conditions.size()>2)
88	{
89	conditions[k-3].ins(0,string_lists[j][i][k]);
90
91	//cout << "modi:" << conditions[k-3] << endl;
92
93	my_rarx->bayes(conditions[k-3]);
94
95
96	//if(k>5)
97	//{
98	// cout << "MaxLik coords:" << emliga->minimal_vertex->get_coordinates() << endl;
99	//}
100
101	}
102
103	}
104
105	//emliga->step_me(0);
106	/*
107	ofstream myfile;
108	myfile.open("c:\\robust_ar1.txt",ios::app);
109	myfile << my_rarx->minimal_vertex->get_coordinates()[0] << ";";
110	myfile.close();
111
112	myfile.open("c:\\robust_ar2.txt",ios::app);
113	myfile << emliga->minimal_vertex->get_coordinates()[1] << ";";
114	myfile.close();
115
116
117	cout << "MaxLik coords:" << emliga->minimal_vertex->get_coordinates() << endl;
118	cout << "Step: " << i << endl;
119	}
120
121	cout << "One experiment finished." << endl;
122
123	ofstream myfile;
124	myfile.open("c:\\robust_ar1.txt",ios::app);
125	myfile << endl;
126	myfile.close();
127
128	myfile.open("c:\\robust_ar2.txt",ios::app);
129	myfile << endl;
130	myfile.close();
131	}*/
132
133
134	// EXPERIMENT: A moving window estimation and prediction of RARX is tested on data generated from
135	// y_t=0.95y_(t-1)+0.05y_(t-2)+0.2*e_t, where e_t is normally, student(4) and cauchy distributed. It
136	// can be compared to the classical setup.
137
138
139	vector<vector<string>> strings;
140
141	char* file_strings[3] = {"c:\\dataADClosePercDiff","c:\\ar_student_single","c:\\ar_cauchy_single"};
142
143	for(int i = 0;i<3;i++)
144	{
145	char dfstring[80];
146	strcpy(dfstring,file_strings[i]);
147	strcat(dfstring,".txt");
148
149	ifstream myfile(dfstring);
150	if (myfile.is_open())
151	{
152	string line;
153	getline(myfile,line);
154
155	vector<string> parsed_line;
156	while(line.find(',') != string::npos)
157	{
158	int loc = line.find(',');
159	parsed_line.push_back(line.substr(0,loc));
160	line.erase(0,loc+1);
161	}
162
163	strings.push_back(parsed_line);
164
165	myfile.close();
166	}
167	}
168
169	for(int j = 0;j<strings.size();j++)
170	{
171	vector<vec> conditions;
172	//emlig* emliga = new emlig(2);
173	RARX* my_rarx = new RARX(3,30,true);
174
175
176	mat V0 = 0.0001 * eye ( 3 );
177	ARX* my_arx = new ARX(0.97);
178	my_arx->set_statistics ( 1, V0 ); //nu is default (set to have finite moments)
179	my_arx->set_constant ( false );
180	my_arx->validate();
181
182
183	for(int k = 1;k<strings[j].size();k++)
184	{
185	vec condition;
186	//condition.ins(0,1);
187	condition.ins(0,strings[j][k]);
188	conditions.push_back(condition);
189
190	//cout << "orig:" << condition << endl;
191
192	if(conditions.size()>1)
193	{
194	conditions[k-2].ins(0,strings[j][k]);
195
196	}
197
198	if(conditions.size()>2)
199	{
200	conditions[k-3].ins(0,strings[j][k]);
201
202	cout << "Condition:" << conditions[k-3] << endl;
203
204	my_rarx->bayes(conditions[k-3]);
205	//my_rarx->posterior->step_me(1);
206
207	vec cond_vec;
208	cond_vec.ins(0,conditions[k-3][0]);
209
210	my_arx->bayes(cond_vec,conditions[k-3].right(2));
211
212
213	if(k>10)
214	{
215	//my_rarx->posterior->step_me(0);
216
217	mat samples = my_rarx->posterior->sample_mat(10);
218
219	pair<vec,mat> imp_samples = my_rarx->posterior->importance_sample(20);
220
221	cout << imp_samples.first << endl;
222
223	vec sample_prediction;
224	for(int t = 0;t<samples.cols();t++)
225	{
226	vec lap_sample = conditions[k-3].left(2);
227	lap_sample.ins(lap_sample.size(),1.0);
228
229	lap_sample.ins(0,LapRNG());
230
231	sample_prediction.ins(0,lap_sample*samples.get_col(t));
232	}
233
234
235	vec sample_pow = sample_prediction;
236
237	// cout << sample_prediction << endl;
238	vec poly_coefs;
239	bool stop_iteration = false;
240	int en = 1;
241	do
242	{
243	double poly_coef = ones(sample_pow.size())*sample_pow/sample_pow.size();
244
245	poly_coef = poly_coefenfact(utility_constant-2+en)/fact(utility_constant-2);
246
247	if(abs(poly_coef)>numeric_limits<double>::epsilon())
248	{
249	sample_pow = elem_mult(sample_pow,sample_prediction);
250	poly_coefs.ins(0,pow(-1.0,en+1)*poly_coef);
251	}
252	else
253	{
254	stop_iteration = true;
255	}
256
257	en++;
258
259	if(en>20)
260	{
261	stop_iteration = true;
262	}
263	}
264	while(!stop_iteration);
265
266	/*
267	ofstream myfile_coef;
268
269	myfile_coef.open("c:\\coefs.txt",ios::app);
270
271	for(int t = 0;t<poly_coefs.size();t++)
272	{
273	myfile_coef << poly_coefs[t] << ",";
274	}
275
276	myfile_coef << endl;
277	myfile_coef.close();
278	*/
279
280	// cout << "Coefficients: " << poly_coefs << endl;
281
282	/*
283	vec bas_coef = vec("1.0 2.0 -8.0");
284	cout << "Coefs: " << bas_coef << endl;
285	cvec actions2 = roots(bas_coef);
286	cout << "Roots: " << actions2 << endl;
287	*/
288
289	cvec actions = roots(poly_coefs);
290
291
292	bool is_max = false;
293	for(int t = 0;t<actions.size();t++)
294	{
295	if(actions[t].imag() == 0)
296	{
297
298
299	double second_derivative = 0;
300	for(int q = 1;q<poly_coefs.size();q++)
301	{
302	second_derivative+=poly_coefs[q]pow(actions[t].real(),q-1)q;
303	}
304
305	if(second_derivative<0)
306	{
307	cout << "Action:" << actions[t].real() << endl;
308
309	is_max = true;
310	}
311	}
312	}
313
314	if(!is_max)
315	{
316	cout << "No maximum." << endl;
317	}
318
319	// cout << "MaxLik coords:" << my_rarx->posterior->minimal_vertex->get_coordinates() << endl;
320
321	double prediction = 0;
322	for(int s = 1;s<samples.rows();s++)
323	{
324
325	double avg_parameter = samples.get_row(s)*ones(samples.cols())/samples.cols();
326
327	prediction += avg_parameter*conditions[k-3][s-1];
328
329
330
331	/*
332	ofstream myfile;
333	char fstring[80];
334	strcpy(fstring,file_strings[j]);
335
336	char es[5];
337	strcat(fstring,itoa(s,es,10));
338
339	strcat(fstring,"_res.txt");
340
341
342	myfile.open(fstring,ios::app);
343
344	//myfile << my_rarx->posterior->minimal_vertex->get_coordinates()[0];
345	myfile << avg_parameter;
346
347	if(k!=strings[j].size()-1)
348	{
349	myfile << ",";
350	}
351	else
352	{
353	myfile << endl;
354	}
355	myfile.close();
356	*/
357	}
358
359	cout << "Prediction: "<< prediction << endl;
360
361	enorm<ldmat>* pred_mat = my_arx->epredictor(conditions[k-3].left(2));
362	double prediction2 = pred_mat->mean()[0];
363
364
365	ofstream myfile;
366	char fstring[80];
367	char f2string[80];
368	strcpy(fstring,file_strings[j]);
369	strcpy(f2string,fstring);
370
371	strcat(fstring,"pred.txt");
372	strcat(f2string,"2pred.txt");
373
374
375	myfile.open(fstring,ios::app);
376
377	// myfile << my_rarx->posterior->minimal_vertex->get_coordinates()[0];
378	myfile << prediction;
379
380	if(k!=strings[j].size()-1)
381	{
382	myfile << ",";
383	}
384	else
385	{
386	myfile << endl;
387	}
388	myfile.close();
389
390
391	myfile.open(f2string,ios::app);
392	myfile << prediction2;
393
394	if(k!=strings[j].size()-1)
395	{
396	myfile << ",";
397	}
398	else
399	{
400	myfile << endl;
401	}
402	myfile.close();
403
404
405	}
406	}
407
408	//emliga->step_me(0);
409	/*
410	ofstream myfile;
411	myfile.open("c:\\robust_ar1.txt",ios::app);
412	myfile << my_rarx->minimal_vertex->get_coordinates()[0] << ";";
413	myfile.close();
414
415	myfile.open("c:\\robust_ar2.txt",ios::app);
416	myfile << emliga->minimal_vertex->get_coordinates()[1] << ";";
417	myfile.close();
418
419
420	cout << "MaxLik coords:" << emliga->minimal_vertex->get_coordinates() << endl;
421	cout << "Step: " << i << endl;*/
422	}
423
424
425	}
426
427
428	// EXPERIMENT: One step ahead price prediction. Comparison of classical and robust model using optimal trading
429	// with maximization of logarithm of one-step ahead wealth.
430
431
432
433	/*
434	cout << "One experiment finished." << endl;
435
436	ofstream myfile;
437	myfile.open("c:\\robust_ar1.txt",ios::app);
438	myfile << endl;
439	myfile.close();
440
441	myfile.open("c:\\robust_ar2.txt",ios::app);
442	myfile << endl;
443	myfile.close();*/
444
445
446	//emlig* emlig1 = new emlig(emlig_size);
447
448	//emlig1->step_me(0);
449	//emlig* emlig2 = new emlig(emlig_size);
450
451	/*
452	emlig1->set_correction_factors(4);
453
454	for(int j = 0;j<emlig1->correction_factors.size();j++)
455	{
456	for(set<my_ivec>::iterator vec_ref = emlig1->correction_factors[j].begin();vec_ref!=emlig1->correction_factors[j].end();vec_ref++)
457	{
458	cout << j << " ";
459
460	for(int i=0;i<(*vec_ref).size();i++)
461	{
462	cout << (*vec_ref)[i];
463	}
464
465	cout << endl;
466	}
467	}*/
468
469	/*
470	vec condition5 = "1.0 1.0 1.01";//"-0.3 1.7 1.5";
471
472	emlig1->add_condition(condition5);
473	//emlig1->step_me(0);
474
475
476	vec condition1a = "-1.0 1.02 0.5";
477	//vec condition1b = "1.0 1.0 1.01";
478	emlig1->add_condition(condition1a);
479	//emlig2->add_condition(condition1b);
480
481	vec condition2a = "-0.3 1.7 1.5";
482	//vec condition2b = "-1.0 1.0 1.0";
483	emlig1->add_condition(condition2a);
484	//emlig2->add_condition(condition2b);
485
486	vec condition3a = "0.5 -1.01 1.0";
487	//vec condition3b = "0.5 -1.01 1.0";
488
489	emlig1->add_condition(condition3a);
490	//emlig2->add_condition(condition3b);
491
492	vec condition4a = "-0.5 -1.0 1.0";
493	//vec condition4b = "-0.5 -1.0 1.0";
494
495	emlig1->add_condition(condition4a);
496	//cout << "************************************************" << endl;
497	//emlig2->add_condition(condition4b);
498	//cout << "************************************************" << endl;
499
500	//cout << emlig1->minimal_vertex->get_coordinates();
501
502	//emlig1->remove_condition(condition3a);
503	//emlig1->step_me(0);
504	//emlig1->remove_condition(condition2a);
505	//emlig1->remove_condition(condition1a);
506	//emlig1->remove_condition(condition5);
507
508
509	//emlig1->step_me(0);
510	//emlig2->step_me(0);
511
512
513	// DA SE POUZIT PRO VYPIS DO SOUBORU
514	// emlig1->step_me(0);
515
516	//emlig1->remove_condition(condition1);
517
518
519
520
521
522	/*
523	for(int i = 0;i<100;i++)
524	{
525	cout << endl << "Step:" << i << endl;
526
527	double condition[emlig_size+1];
528
529	for(int k = 0;k<=emlig_size;k++)
530	{
531	condition[k] = (rand()-RAND_MAX/2)/1000.0;
532	}
533
534
535	vec* condition_vec = new vec(condition,emlig_size+1);
536	emlig1->add_condition(*condition_vec);
537
538	/*
539	for(polyhedron* toprow_ref = emlig1->statistic.rows[emlig_size]; toprow_ref != emlig1->statistic.end_poly; toprow_ref = toprow_ref->next_poly)
540	{
541	cout << ((toprow*)toprow_ref)->probability << endl;
542	}
543	*/
544	/*
545	cout << emlig1->statistic_rowsize(emlig_size) << endl << endl;
546
547	/*
548	if(i-emlig1->number_of_parameters >= 0)
549	{
550	pause(30);
551	}
552	*/
553
554	// emlig1->step_me(i);
555
556	/*
557	vector<int> sizevector;
558	for(int s = 0;s<=emlig1->number_of_parameters;s++)
559	{
560	sizevector.push_back(emlig1->statistic_rowsize(s));
561	}
562	*/
563	//}
564
565
566
567
568	/*
569	emlig1->step_me(1);
570
571	vec condition = "2.0 0.0 1.0";
572
573	emlig1->add_condition(condition);
574
575	vector<int> sizevector;
576	for(int s = 0;s<=emlig1->number_of_parameters;s++)
577	{
578	sizevector.push_back(emlig1->statistic_rowsize(s));
579	}
580
581	emlig1->step_me(2);
582
583	condition = "2.0 1.0 0.0";
584
585	emlig1->add_condition(condition);
586
587	sizevector.clear();
588	for(int s = 0;s<=emlig1->number_of_parameters;s++)
589	{
590	sizevector.push_back(emlig1->statistic_rowsize(s));
591	}
592	*/
593
594	return 0;
595	}
596
597

Note: See TracBrowser for help on using the browser.

Context Navigation

root/applications/robust/main.cpp @ 1343

Download in other formats: