solutions for fitting exercises

statistics/code/boltzmann.m

@@ -0,0 +1,7 @@
+function y = boltzmann(parameter, x)
+% parameter 1: alpha
+% parameter 2: k
+% parameter 3: x_0
+% parameter 4: y_0
+
+y = (parameter(1) ./ (1 + exp(-parameter(2) .* (x - parameter(3))))) + parameter(4);

statistics/code/create_linear_data.m

@@ -0,0 +1,7 @@
+function y = create_linear_data(x)
+
+m = 2.5;
+n = -0.35;
+d = 2.5;
+
+y = x .* m + n + randn(size(x)) .* d;

statistics/code/estimate_regression.m

@@ -0,0 +1,7 @@
+function param = estimate_regression(x,y,p_0)
+
+objective_function = @(p)lsq_error(p, x, y);
+param = fminunc(objective_function, p_0);
+disp(param)
+param1 = fminsearch(objective_function, p_0);
+disp(param1)

statistics/code/exponential.m

@@ -0,0 +1,5 @@
+function y = exponential(parameter, x)
+% Function implements an exponential function with two parameters
+% controlling the amplitude and the time constant.
+
+y = parameter(1) .* exp(x./parameter(2));

statistics/code/lsq_gradient_sigmoid.m

@@ -0,0 +1,9 @@
+function gradient = lsq_gradient_sigmoid(parameter, x, y)
+h = 1e-6;
+
+gradient = zeros(size(parameter));
+for i = 1:length(parameter)
+    parameter_h = parameter;
+    parameter_h(i) = parameter_h(i) + h;
+    gradient(i) = (lsq_sigmoid_error(parameter_h, x, y) - lsq_sigmoid_error(parameter, x, y)) / h;
+end

statistics/code/lsq_sigmoid_error.m

@@ -0,0 +1,8 @@
+function error = lsq_sigmoid_error(parameter, x, y)
+% p(1) the amplitude
+% p(2) the slope
+% p(3) the x-shift
+% p(4) the y-shift
+
+y_est = parameter(1)./(1+ exp(-parameter(2) .* (x - parameter(3)))) + parameter(4);
+error = mean((y_est - y).^2);

statistics/code/sigmoidal_gradient_descent.m

@@ -0,0 +1,44 @@
+
+
+%% fit the sigmoid
+
+clear 
+close all
+
+load('iv_curve.mat')
+
+figure()
+plot(voltage, current, 'o')
+xlabel('voltate [mV]')
+ylabel('current [pA]')
+
+% amplitude, slope, x-shift, y-shift
+%parameter = [10 0.25 -50, 2.5];
+parameter = [20 0.5 -50, 2.5];
+
+eps = 0.1;
+% do the descent
+gradient = [];
+steps = 0;
+error = [];
+
+while isempty(gradient) || norm(gradient) > 0.01
+    steps = steps + 1; 
+    gradient = lsq_gradient_sigmoid(parameter, voltage, current);
+    error(steps) = lsq_sigmoid_error(parameter, voltage, current);
+    parameter = parameter - eps .* gradient; 
+end
+plot(1:steps, error)
+
+disp('gradient descent done!')
+disp(strcat('final position: ', num2str(parameter)))
+disp(strcat('final error: ', num2str(error(end))))
+
+%% use fminsearch
+parameter = [10 0.5 -50, 2.5];
+
+objective_function = @(p)lsq_sigmoid_error(p, voltage, current);
+param = fminunc(objective_function, parameter);
+disp(param)
+param1 = fminsearch(objective_function, parameter);
+disp(param1)