diff --git a/projects/project_input_resistance/solution/exponentialDecay.m b/projects/project_input_resistance/solution/exponentialDecay.m new file mode 100644 index 0000000..2ee27cf --- /dev/null +++ b/projects/project_input_resistance/solution/exponentialDecay.m @@ -0,0 +1,5 @@ +function y = exponentialDecay(p, x) + + y = p(1) .* exp(-x./p(2)); + + diff --git a/projects/project_input_resistance/solution/gradientDescent.m b/projects/project_input_resistance/solution/gradientDescent.m new file mode 100644 index 0000000..867e863 --- /dev/null +++ b/projects/project_input_resistance/solution/gradientDescent.m @@ -0,0 +1,14 @@ +function [position, errors] = gradientDescent(x, y, position) + +gradient = []; +errors = []; +count = 1; +eps = 0.02; +close all +while isempty(gradient) || norm(gradient) > 0.025 + gradient = lsqGradient(position, x,y); + % disp(gradient) + errors(count) = lsqError(position, x, y); + position = position - eps .* gradient'; + count = count + 1; +end \ No newline at end of file diff --git a/projects/project_input_resistance/solution/lsqError.m b/projects/project_input_resistance/solution/lsqError.m new file mode 100644 index 0000000..1a48a52 --- /dev/null +++ b/projects/project_input_resistance/solution/lsqError.m @@ -0,0 +1,12 @@ +function error = lsqError(parameter, x, y) +% Objective function for fitting a exponential equation to data. +% +% Arguments: parameter, vector containing the parameters for the exp. decay +% x, vector of the input values +% y, vector of the corresponding measured output values +% +% Returns: the estimation error in terms of the mean sqaure error + + y_est = exponentialDecay(parameter, x); + error = meanSquareError(y, y_est); +end diff --git a/projects/project_input_resistance/solution/lsqGradient.m b/projects/project_input_resistance/solution/lsqGradient.m new file mode 100644 index 0000000..44b3881 --- /dev/null +++ b/projects/project_input_resistance/solution/lsqGradient.m @@ -0,0 +1,9 @@ +function gradient = lsqGradient(parameter, x, y) + + h = 1e-1; % stepsize for derivatives + gradient = zeros(length(parameter),1); + for i = 1:length(parameter) + new_p = parameter; + new_p(i) = parameter(i) + h; + gradient(i) = (lsqError(new_p, x, y) - lsqError(parameter, x, y))/ h; + end diff --git a/projects/project_input_resistance/solution/meanSquareError.m b/projects/project_input_resistance/solution/meanSquareError.m new file mode 100644 index 0000000..253d31a --- /dev/null +++ b/projects/project_input_resistance/solution/meanSquareError.m @@ -0,0 +1,10 @@ +function error = meanSquareError(y, y_est) +% Mean squared error between observed and predicted values. +% +% Arguments: y, vector of observed values. +% y_est, vector of predicted values. +% +% Returns: the error in the mean-squared-deviation sense. + + error = mean((y - y_est).^2); +end diff --git a/projects/project_input_resistance/solution/project_input_resistance.m b/projects/project_input_resistance/solution/project_input_resistance.m new file mode 100644 index 0000000..5e52441 --- /dev/null +++ b/projects/project_input_resistance/solution/project_input_resistance.m @@ -0,0 +1,43 @@ +clear +maps = {'CMS', 'MS'}; + +taus = zeros(length(dir('data/resistance/resistance*.mat')), 1); +cell_types = cell(size(taus)); +count = 1; +size(taus) +for j = 1:length(maps) + files = dir(strcat('data/resistance/resistance_', maps{j}, '*.mat')); + for i = 1:length(files) + load(strcat('data/resistance/', files(i).name)); + disp(files(i).name) + x = t(t>0.5 & t<15); + x = x - x(1); + y = V(t>0.5 & t<15); + y = y - y(end); + initial_params = [1, 2.5]; + [params, errors] = gradientDescent(x, y, initial_params); + disp(count) + taus(count) = params(2); + cell_types{count} = maps{j}; + count = count + 1; + plot_fit(x, y, errors, params); + end +end + + +function plot_fit(x, y, errors, parameter) + figure() + subplot(2,1,1) + hold on + scatter(x,y, 'displayname', 'data') + f_x = exponentialDecay(parameter, x); + plot(x, f_x, 'displayname', 'fit') + xlabel('Input') + ylabel('Output') + grid on + legend show + subplot(2,1,2) + plot(errors) + xlabel('optimization steps') + ylabel('error') +end