MATH590_C1.3

MATH590 Homework 1: Approximation in R d Carry out analysis of the cell-phone accelerometer data according to the following template, answering the questions. The objective of this homework is to familiarize yourself with some linear algebra and numerical analysis techniques that are commonly encountered in data analysis. Turn in this TeXmacs file and the two data files AveragePeriod.LastName.dat a, PolyCoef.LastName.data . The results contained in these ±les will be used in subsequent topological data analysis. 1 Qualitative data analysis A ±rst step in processing the data ( a ; " ) acquired from the cell phone is to carry out some basic qualitative analysis, shown here using Octave (Matlab clone). 1.1 Data input 1.1.1 Visual inspection of data octave> dir='/home/student/courses/MATH590/NUMdata'; chdir(dir); LastName = 'Mitran'; data=csvread(strcat(LastName,'.csv')); [m,d] = size(data); disp([m d]); 1

18841 8 octave> data(1:30,:) 0 B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B B @ 0 0 0 0 0 0 0 0 0.004 ¡ 0.7803 0.2999 ¡ 1.799 ¡ 0.1402 0.2944 0.0533 0 0.005 ¡ 0.7803 0.2999 ¡ 1.799 0.0015 0.3671 0.2121 0 0.005 ¡ 0.0953 0.0762 ¡ 0.0111 0.0015 0.3671 0.2121 0 0.006 ¡ 0.0953 0.0762 ¡ 0.0111 0.0015 0.3671 0.2121 0 0.006 ¡ 0.0953 0.0762 ¡ 0.0111 0.0015 0.3671 0.2121 0 0.006 ¡ 0.0953 0.0762 ¡ 0.0111 0.0015 0.3671 0.2121 0 0.006 ¡ 0.0953 0.0762 ¡ 0.0111 0.0015 0.3671 0.2121 0 0.006 ¡ 0.0953 0.0762 ¡ 0.0111 0.0015 0.3671 0.2121 0 0.007 ¡ 0.0953 0.0762 ¡ 0.0111 ¡ 0.0003 0.7379 ¡ 0.0279 0 0.007 ¡ 0.6242 ¡ 0.0203 ¡ 0.9742 ¡ 0.0003 0.7379 ¡ 0.0279 0 0.025 ¡ 0.6242 ¡ 0.0203 ¡ 0.9742 ¡ 0.0003 0.7379 ¡ 0.0279 0 0.056 ¡ 0.6242 ¡ 0.0203 ¡ 0.9742 ¡ 0.0003 0.7379 ¡ 0.0279 0 0.056 ¡ 0.6242 ¡ 0.0203 ¡ 0.9742 ¡ 0.0003 0.7379 ¡ 0.0279 0 0.056 ¡ 0.6242 ¡ 0.0203 ¡ 0.9742 ¡ 0.0003 0.7379 ¡ 0.0279 0 0.068 ¡ 0.6242 ¡ 0.0203 ¡ 0.9742 ¡ 0.0003 0.7379 ¡ 0.0279 0 0.068 ¡ 0.6242 ¡ 0.0203 ¡ 0.9742 ¡ 0.1671 0.3402 0.2799 0 0.068 ¡ 1.0653 0.5028 ¡ 0.6942 ¡ 0.1671 0.3402 0.2799 0 0.101 ¡ 1.0653 0.5028 ¡ 0.6942 ¡ 0.1671 0.3402 0.2799 0 0.119 ¡ 1.0653 0.5028 ¡ 0.6942 ¡ 0.1671 0.3402 0.2799 0 0.12 ¡ 1.0653 0.5028 ¡ 0.6942 ¡ 0.1671 0.3402 0.2799 0 0.12 ¡ 1.0653 0.5028 ¡ 0.6942 ¡ 0.1671 0.3402 0.2799 0 0.126 ¡ 1.0653 0.5028 ¡ 0.6942 ¡ 0.1671 0.3402 0.2799 0 0.126 ¡ 1.0653 0.5028 ¡ 0.6942 ¡ 0.0082 0.5577 0.2268 0 0.127 ¡ 0.1956 0.0971 ¡ 0.7178 ¡ 0.0082 0.5577 0.2268 0 0.176 ¡ 0.1956 0.0971 ¡ 0.7178 ¡ 0.0082 0.5577 0.2268 0 0.177 ¡ 0.1956 0.0971 ¡ 0.7178 ¡ 0.0082 0.5577 0.2268 0 0.177 ¡ 0.1956 0.0971 ¡ 0.7178 ¡ 0.0082 0.5577 0.2268 0 0.186 ¡ 0.1956 0.0971 ¡ 0.7178 ¡ 0.0082 0.5577 0.2268 0 0.186 ¡ 0.1956 0.0971 ¡ 0.7178 ¡ 0.1506 0.5546 0.3117 0 1 C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C C A octave> 2

gait. The largest angular velocities correspond to yaw, also as expected since this is the alternating forward-backward motion of left-right side of the body during the gait. Finally, walking is an e±ort against gravitational acceleration g = 9.8 m / s 2 , and one would expect values about one-tenth of g . 1.1.2 Plot data Write data to a ±le. octave> [t,iu]=unique(data(:,1)); mu=max(size(iu)); raddeg=pi/180; a=data(iu,2:4); epsilon=data(iu,5:7)/raddeg; fid=fopen(strcat(LastName,'.data'),'w'); fprintf(fid,'%f %f %f %f %f %f %f\n',[t a epsilon]'); fclose(fid); octave> Change the ±le plotted to your data. GNUplot] cd '/home/student/courses/MATH590/NUMdata' set terminal postscript eps enhanced color set style line 1 lt 2 lc rgb "red" lw 3 set style line 2 lt 2 lc rgb "green" lw 3 set style line 3 lt 2 lc rgb "blue" lw 3 plot 'Mitran.data' u 1:2 w l ls 1 title "a1", '' u 1:3 w l ls 2 title "a2", '' u 1:4 w l ls 3 title "a3" -4 -3 -2 -1 0 1 2 3 4 5 0 20 40 60 80 100 120 140 160 180 a1 a2 a3 4

GNUplot] GNUplot] cd '/home/student/courses/MATH590/NUMdata' set terminal postscript eps enhanced color set style line 1 lt 2 lc rgb "red" lw 3 set style line 2 lt 2 lc rgb "green" lw 3 set style line 3 lt 2 lc rgb "blue" lw 3 plot 'Mitran.data' u 1:5 w l ls 1 title "epsilon1", '' u 1:6 w l ls 2 title "epsilon2", '' u 1:7 w l ls 3 title "epsilon3" -150 -100 -50 0 50 100 150 0 20 40 60 80 100 120 140 160 180 epsilon1 epsilon2 epsilon3 GNUplot] Question 2. What are your observations on the physical relevance of the data? Data shows an average close to zero, as expected, with large values for yaw angular velocity 5

octave> GNUplot] cd '/home/student/courses/MATH590/NUMdata' set terminal postscript eps enhanced color set style line 1 lt 2 lc rgb "green" lw 3 plot "InterpolatedVerticalAcceleration.data" u 1:2 w l ls 1 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 2.5 3.5 4 4.5 5 5.5 6 6.5 "InterpolatedVerticalAcceleration.data" u 1:2 GNUplot] Question 3. How do you explain changes in the vertical acceleration between steps? This is the expected acceleration of the body's center of mass during stepping motion. In contrast if this was rolling motion (i.e., on wheels) the vertical acceleration should be small. This values shows even large variation when stepping on stairs. 2.2 Determine gait period octave> AGi=fft(aGi); PAGi = log10(AGi.*conj(AGi)); fid=fopen('GaitPowerSpectrum.data','w'); fprintf(fid,'%f\n',PAGi(1:nG/2-1)); fclose(fid); octave> GNUplot] cd '/home/student/courses/MATH590/NUMdata' set terminal postscript eps enhanced color set style line 1 lt 2 lc rgb "green" lw 3 plot 'GaitPowerSpectrum.data' w l ls 1 7

-2 -1 0 1 2 3 4 5 0 20 40 60 80 100 120 140 ’GaitPowerSpectrum.data’ GNUplot] octave> [val,idx] = max(PAGi); disp([val idx]); 4.6046 6.0000 octave> TG = (tG1-tG0)/(idx-1) 0.5862 octave> nT=floor(max(size(aGi))/(idx-1)) 51 octave> Question 4. Does the period value correspond to a full stride or stepping on one leg? The period corresponds to stepping on one leg. A full stride (returning to the same foot on the ground at start of step) corresponds to 2 periods. 2.3 Splice data into multiple periods Find peak vertical accelerations. octave> [aPeak, iPeak] = findpeaks(aGi-min(aGi),"MinPeakHeight",1., "MinPeakDistance",nT/3); fid=fopen('TrajectoryPeaks.data','w'); ta = [tGi(iPeak)' (aPeak+min(aGi))']; fprintf(fid,'%f %f\n',ta'); fclose(fid); octave> nPeriods = max(size(iPeak))-1 8

0.53811 0.59536 0.60681 0.53811 octave> i=1; while(i<=nPeriods) tT(1:nT(i),i) = tT(1:nT(i),i) - tT(1,i); tT(1:nT(i),i) = tT(1:nT(i),i)/TT(i); amp = max(aT(1:nT(i),i)) - min(aT(1:nT(i),i)); aT(1:nT(i),i) = aT(1:nT(i),i)/amp; i++; endwhile; octave> nTi=100; dt=1./(nTi-1); ti=(0:nTi-1)*dt; ti=ti'; aTi=zeros([nTi,nPeriods]); g=zeros([nTi,1]); i=1; while(i<=nPeriods) ai = interp1(tT(1:nT(i),i),aT(1:nT(i),i),ti',"nearest"); aTi(:,i) = ai; g = g + ai'; i++; endwhile; g = g/nPeriods; g = g/(max(g)-min(g)); octave> fid=fopen('Periods.data','w'); ta = [ti aTi g]; fprintf(fid,'%f %f %f %f %f %f\n',ta'); fclose(fid); fid=fopen(strcat(strcat('AveragePeriod.',LastName),'.data'),'w'); ta = [ti g]; fprintf(fid,'%f %f \n',ta'); fclose(fid); octave> GNUplot] cd '/home/student/courses/MATH590/NUMdata' set terminal postscript eps enhanced color set style line 1 lt 2 lc rgb "red" lw 6 set style line 2 lt 2 lc rgb "green" lw 3 set style line 3 lt 2 lc rgb "blue" lw 3 plot "Periods.data" u 1:2 w l ls 3, '' u 1:3 w l ls 3, '' u 1:4 w l ls 3, '' u 1:5 w l ls 3, '' u 1:6 w l ls 1 10

-0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 0 0.2 0.4 0.6 0.8 1 "Periods.data" u 1:2 ’’ u 1:3 ’’ u 1:4 ’’ u 1:5 ’’ u 1:6 GNUplot] Question 5. What are the potential drawbacks of de±ning an ²average³ gait? Consider the limiting cases of too few or very many sample periods. For too many periods, di´erences in the path (climbing, turning, descending) can mask the average gait assumed to be speci±c to a person. For, say, a single period, the gait might not be typical, e.g., due to a cough or irregularity in the path. It would be best to more carefully control the path in order to identify the person, e.g., limit data to walking a straight path of 10 meters. 2.4 Least squares Seek a more economical representation of the average gait waveform, currently stored as a vector g 2 R m of the vertical acceleration values at times within the vector t 2 R m . For example, consider approximating the waveform by a parabola p ( t ) = c 0 + c 1 t + c 2 t 2 , leading to the least squares problem min c 2 R 3 k Lc ¡ g k ; L = ¡ 1 t t 2 ± 2 R m ± 3 ; t k = ¡ t 1 k ::: t m k ± T : (1) A solution is found by projection onto the column space of L , QR = L ; P C ( L ) = QQ T ; P C ( L ) g = QRc ) Rc = Q T g : (2) octave> L=[ti.^0 ti ti.^2]; [Q,R]=qr(L,0); [size(Q); size(R)] 11

2.5 Min-max An alternative representation is through a linear combination of Chebyshev polynomials, q ( t ) = d 0 T 0 ( t ) + d 1 T 1 ( t ) + d 2 T 2 ( t ) known to be a good approximation of the min-max polynomial of the data. The coe±cient vector d 2 R 3 is more di±cult to ²nd by comparison to the least squares case, and is carried out through a procedure known as the exchange algorithm, implemented in Octave by the polyfitinf function. octave> d=polyfitinf(2,nTi,0,ti,g,5.0E-7,nTi) ( 3.3996 ¡ 3.665 0.74003 ) octave> q = polyval(d,ti); fid=fopen('MinMax.data','w'); ta = [ti g q]; fprintf(fid,'%f %f %f\n',ta'); fclose(fid); octave> GNUplot] cd '/home/student/courses/MATH590/NUMdata' set terminal postscript eps enhanced color set style line 1 lt 2 lc rgb "red" lw 6 set style line 2 lt 2 lc rgb "green" lw 3 set style line 3 lt 2 lc rgb "blue" lw 6 plot "LeastSquares.data" u 1:3 w l ls 1, '' u 1:2 w p ls 2, "MinMax.data" u 1:3 w l ls 3 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 "LeastSquares.data" u 1:3 ’’ u 1:2 "MinMax.data" u 1:3 GNUplot] 13

Question 7. Which gait approximation is better suited to walker identi±cation, least squares or min-max? Explain your reasoning. Min-max is suited to approximation that minimizes a presumed ²maximal³ error. The problem is that such a maximal error cannot be de±ned in this case, whereas, for example, it is straightforward to de±ne in the approximation of an analytical function such as sin(x). The least-squares approximation is better in this case. Save the coe±cients of the two representations (least-squares and Chebyshev) octave> fid=fopen(strcat(strcat('PolyCoef.',LastName),'.data'),'w'); ta = [c d']; fprintf(fid,'%f %f\n',ta'); fclose(fid); octave> 14