Recently, Steve wrote a blog discussing code clarity/obscurity in the context of one-line code solutions. Simply stated, the problem he solved is this. Find the largest value in an array adjacent to a zero value.
Must be zillions of ways, right?
Let me start with an example. If the array in question is A,
A = [1, 5, 3, 0, 2, 7, 0, 8, 9, 1 0];
the correct answer is 8.
Steve's First Solution
Steve, being an image processing guru, immediately arrived at a solution using image processing techniques. For me, not so
obvious what is going on! If I inherited this code, I'd have some work cut out for me to understand it.
Steve then noticed his solution was only good for positive values.
B = [5 4 -1 0 -2 0 -5 8];
fs(B)
ans =
0
Steve's Second Answer
So Steve tried again. Again, though I understand what a mask is, I don't fully comprehend the solution, and would probably
resort to reading the imdilate documentation from Image Processing Toolbox.
Here we see Doug's answer, written without dependence on anything but MATLAB. Very nice solution, but I have to think about
it a bit. Again, I know what a mask is, I see we are looking for NaN values after division by 0, so we'll find the zeros, convolve to get the neighbors of the zeros, and find the maximum of
those. Whew! I'm working hard now but I get it.
The next suggestion that appeared is from Jos. It's dense, but ultimately comprehensible to me. First, create a 2 row vector
with the data padded at each end with the value 1, and aligned so the first value is above the 3rd, the 2nd above the 4th,
etc. (i.e., shifted over by 2 elements). Check each column to see if any of the values are 0. Select the data from the columns
with a zero value, and find the maximum of those values.
fj =
@(X)max(X(any([1,X(1:end-1);X(2:end),1]==0)))
ans =
8
ans =
-1
Start Over
My turn to start over and think about the problem. First, I want to find the indices of zeros in the array. From there,
I construct an array of the neighbors.
zind = find(B==0)
maxcandidates = [zind-1 zind+1]
zind =
4 6
maxcandidates =
3 5 5 7
Next, I weed out values outside the range (possibly the first and last).
Now I index into my candidates and get the maximum value!
max(B(maxcandidates))
ans =
-1
As Chef Tell would say, "Very simple, very easy."
I could write this algorithm in one statement, but I fear it would immediately lose clarity.
Repeated Zeros - Clarity of Problem Statement
I got to wondering what result the originator of the question had in mind if there were repeated zeros. Why? Because I could
possibly get a zero result? Was that intended? All the solutions here allow that, so mine does too and I don't need to do
anything more. Let's check what happens.
C = [1 2 -4 0 0];
fj(C)
ans =
0
Your Thoughts?
It can definitely be fun and challenging to write nifty, compact, and very dense one-liners in MATLAB. But the horror of
going back to that code later, deciphering it so it can be updated or whatever, is not for me typically. I like to write
code where the comments don't have to be a lot longer than the code! What are your thoughts? Post them here.
I recently wrote a pair of posts (1 and 2) about finding roots of equations, and how you might use MATLAB to explore this topic with emphasis on the method of fixed point iteration.
Let me restate the example function. Let's start with a simple cubic polynomial f(x)
which we define using an anonymous function.
f = @(x) x.^3 + x - 1;
fplot(f,[-2 2])
title f
grid on
First Fixed Point Iteration Attempt
Previously we rewrote f(x)=0 so that x was alone on one side of the equation.
g1 = @(x) 1 - x.^3;
This function, g1(x) did not help find the root between 0 and 1 - every step took us further away from the solutions we found with roots and fzero.
fzsolution = fzero(f,0.5)
fzsolution =
0.68233
Second Way to Rewrite Equation
There's another way we can write the equation for the fixed point. Remember, we want to rearrange f(x)=0 into something like g2(x)=x. We have already tried g1 = 1 - x^3. We can also rewrite the equation as
g2 = @(x) 1./(x.^2+1);
fplot(g2,[0 1]);
hold on
straightLine = @(x) x;
fplot(straightLine, [0 1], 'g')
legend('g2','x','Location','SouthEast')
grid on
axis equal, axis([0 1 0 1])
Let's Iterate
Following the same idea from the last post, we now iterate, finding successive guesses of x, computing the value of g2(x), setting this value to be the next estimate of x, etc.
Let's Try to "Zero" in on the Solution
First, we'll select a starting value x = 0.5, and y, from the straight line also at 0.5.
And let's plot the trajectory of the solutions, starting with the first point.
plot(x,y,'r',x(1),y(1),'r*')
First Iteration
As described in the previous post, I now set the current value of g2 to the new x value, sliding onto the straight line, and then calculate the next value of g2 from this new value for x. And plot it.
Let's do a few more iterations and preallocate the arrays instead of growing them. And yes, I know in this post I kept overplotting
lines successively. I didn't want us to get distracted by managing the plots.
x(9:22) = NaN;
y(9:22) = NaN;
for n = 9:2:21
x(n) = y(n-1);
y(n) = y(n-1);
x(n+1) = x(n);
y(n+1) = g2(x(n+1));
end
plot(x,y,'r')
hold off
Fixed Point Found!
We appear to be converging on a fixed point, since after iterating, the final values for x and g2(x) (which is y) are getting quite close to each other.
final = [ x(end) y(end) fzsolution]
final =
0.68039 0.68356 0.68233
Wrapping Up the Series of Posts
This post completes this series of posts on fixed point iteration. There is, of course, more that you could do in a class.
For example, you might explore what characteristics of the functions g1 and g2 made the fixed point iteration strategy behave differently. Perhaps look at first derivatives?
Is this sort of tutorial relevant to your work, especially if you are teaching? How about the incorporation of the visualization
during the exploration? Let me know your thoughts here.
In a recent post, I started a series of posts about finding roots of equations, and how you might use MATLAB to teach some of this topic.
In particular, I will introduce the idea of fixed point iteration.
Let me restate the example function. We'll start with a simple cubic polynomial f
which I define as an anonymous function.
f = @(x) x.^3 + x - 1;
fplot(f,[-2 2])
title f
grid on
Fixed Point Algorithm
One way to find a solution to f(x) = 0 is to define another function, g such that g(x) = x where the function g is related to f. If we can find a value x that solves the equation with g, this fixed point is also the zero of the function f. How do we define such a suitable g? Looking at the equation for f and setting its value to zero, we can derive an equation for x.
Here's one way we can rewrite f(x) = 0 so that x is alone on one side of the equation.
Let's call this function g1.
g1 = @(x) 1 - x.^3;
fplot(g1,[-2 2]);
g1str = '1-x^3';
title g1
grid on
Now let's add a straight line to represent
straightLine = @(x) x;
title ''
hold on
fplot(straightLine, [-2 2], 'g')
legend(g1str,'x','Location','SouthEast')
axis([-2 2 -2 2]), axis equal
The intersection of these 2 curves is a fixed point of g1(x), and therefore a zero of f(x). Let's zoom in a bit.
axis([0 1 0 1])
That intersection looks like it happens around the same location as the answers we got from using roots and fzero in the previous post.
fzsolution = fzero(f,0.5)
fzsolution =
0.68233
plot(fzsolution, g1(fzsolution),'b*')
Let's Try to "Zero" in on the Solution
First, we'll select the same starting value that we used when we called fzero. The starting point is then (0.5,0.5). And then we can easily calculate the next value of g1.
Let's plot the trajectory of the solutions, noting the first point with a red star.
plot(x,y,'r',x(1),y(1),'r*')
First Iteration
Now let's get the next estimate. We'll set the current value of g1 to the new x value, sliding onto the straight line, and then calculate the next value of g1 from this new value for x. And plot it.
n = 7;
x(n) = y(n-1);
y(n) = y(n-1);
x(n+1) = x(n);
y(n+1) = g1(x(n+1));
plot(x,y,'r')
hold off
The Conundrum
By now, you have surely noticed that instead of spiraling into the solution we found with fzero, our estimates are spiraling out! In the next post, we'll explore another definition for g(x) = x and see what happens from there.
Confession
Yes, I know in this post I kept overplotting lines successively. But I didn't want us to get distracted by managing the plots.
The Series of Posts
In addition to this post and the previous one, there will be one more. Until it is published, the 3rd link will not be useful.
I've been interested in teaching for a long time, including ways to use MATLAB. One concept that students might need to understand
early in their college careers is that of finding roots or zeros of functions. To understand at least some of the algorithms,
you might want to teach the students about fixed points for functions. It's the basis for some methods of solving equations or finding roots, algorithms such as Newton's method, finding square roots, and more.
Have you ever tried to overlap more than one plot, only to find that you had to fiddle with a few times before it looked right?
Perhaps this is because you aren't take full advantage of the function hold.
Suppose I want to plot two sets of data or functions on the same plot. How might I do that? There are several ways. Let's
create some data and a function.
t = 0:0.005:1;
f = @(t) sin(2*pi*10*t);
s = f(t);
sn = s + 0.3*randn(size(t));
Plot the Data and the Function
If I have everything assembled up front, I can plot all the values at once. The function plot is smart enough to cycle through colors for multiple lines.
plot(t,s,t,sn)
legend('signal','signal with noise','Location','SouthEast')
Try Another Way
Suppose we have the first data values.
plot(t,s)
And then need to add some more values later. Here's a way I can do this.
hold on
plot(t,sn)
legend('signal','signal with noise','Location','SouthEast')
hold off
But the plot is a bit confusing, no? The legend isn't very helpful since each plot has the same color and linestyle. Instead, I can specify the color(s) I want with the
plot commands.
plot(t,s)
hold on
plot(t,sn,'g')
legend('signal','signal with noise','Location','SouthEast')
hold off
That green is too bright. To use the one from the first plot, I can find out the correct green color by looking at the defaults.
You can see that the green is definitely toned down (the color represented by RGB values in the allcolors(2,:).
plot(t,s)
hold on
plot(t,sn,'color',allcolors(2,:))
legend('signal','signal with noise','Location','SouthEast')
hold off
Why Use Hold?
Why use hold at all if it's giving me headaches? Because sometimes I want to add data to the plot using a different plot function. Here's
an example.
fplot(f, [0 1])
hold on
plot(t,sn,'color',allcolors(2,:))
legend('signal','signal with noise','Location','SouthEast')
hold off
A Better Way
Surely there must be a better way than managing the colors on your own. And there is! The all option to hold magically (quoting from the doc)
holds the plot and the current line color and line style so that subsequent plotting commands do not reset
the ColorOrder and ColorOrder property values to the beginning of the list. Plotting commands continue cycling through the
predefined colors and linestyles from where the last plot stopped in the list.
Let's try it.
fplot(f, [0 1])
hold all
plot(t,sn)
legend('signal','signal with noise','Location','SouthEast')
hold off
NOTE: I updated this since I didn't use hold off in my examples earlier so old lines were overwriting new ones in the legend. Much better now.
Will hold all Help You?
I'd be curious to hear if the behavior of hold all helps you out. Let me know here.
When I got to work last Friday, I saw an email discussion, on behalf of a customer, trying to find a good way to add a new
field to a struct array. So this post will start with that problem, and then show a different way to collect the same information, in a dataset array.
Let's quickly check that we get the same results with each technique.
allsame = isequal(S1,S2,S3,S4)
allsame =
1
What's the Data Look Like?
It's hard to look at the data here (in, e.g., S1) because the contents of each struct element is completely at the users's disposal. So I can look at one array element at a time.
S1(1)
ans =
Name: 'John'
Age: 26
Height: 168
Or I can look at all of the data in a single field at once.
[S1.Age]
ans =
26 18
But I don't get to see all of the data in one glance.
Completely Different View
And now for something completely different. I've blogged before about dataset arrays from Statistics Toolbox. Here's another instance where one might be useful. I treat the columns like individual fields, and the rows as individual
records. Each column contains data of a single datatype. Here's the data.
names =
'John'
'Henri'
d1 =
Name Age
'John' 26
'Henri' 18
Two things to note here in contrast to using a struct to contain the information. First, the arguments appear in a different order in the two solutions. Second, the numeric data
doesn't need to be placed in a cell array for the dataset, making the data management more natural, in my opinion.
Let me make a new dataset with additional data, heights.
Now let me collect the original dataset d1 with the new information in d2. Here are some ways to achieve this. First, just use square brackets ([]) as you would for regular array concatenation.
dnew1 = [d1 d2]
dnew1 =
Name Age Height
'John' 26 168
'Henri' 18 175
Another way to do this is to add the information in a struct-like way to the original dataset.
dnew2 = d1;
dnew2.Height = [168; 175]
dnew2 =
Name Age Height
'John' 26 168
'Henri' 18 175
Now let's make different dataset with new information, but with the order of the 2 entries swapped.
What happens if we try to collect d1 and d3 together into one dataset?
try
dnew3 = [d1 d3];
catch ExcDataset
disp(ExcDataset.message)
end
Duplicate variable names with distinct data.
As you can see, I can't just collect them together via concatenation. However, I can combine or join the two datasets correctly.
dnew3 = join(d1,d2,'Name')
dnew3 =
Name Age Height
'John' 26 168
'Henri' 18 175
Notice how easily I can see all the data at once here, compared to the struct array.
How Do You Arrange Your Data?
Do you use either of these strategies for arranging your data (struct or dataset arrays)? Or do you do something different? I'd love to hear your experiences here.
With a little more thought, I was able to adapt the last program I showed to incorporate empty values taking on the defaults.
This allows the user to override defaults later in the argument list without having to reset all the previous values. Let's
see how this works.
type somefun2AltEmptyDefs
function y = somefun2AltEmptyDefs(a,b,varargin)
% Some function that requires 2 inputs and has some optional inputs.
% only want 3 optional inputs at most
numvarargs = length(varargin);
if numvarargs > 3
error('myfuns:somefun2Alt:TooManyInputs', ...
'requires at most 3 optional inputs');
end
% set defaults for optional inputs
optargs = {eps 17 @magic};
% skip any new inputs if they are empty
newVals = cellfun(@(x) ~isempty(x), varargin);
% now put these defaults into the valuesToUse cell array,
% and overwrite the ones specified in varargin.
optargs(newVals) = varargin(newVals);
% or ...
% [optargs{1:numvarargs}] = varargin{:};
% Place optional args in memorable variable names
[tol, mynum, func] = optargs{:};
Compare Functions
I can use the Compare Against tool (from the Tools menu) to compare my new function with my previous one. Here's a screenshot.
You can click on it to see a legible version.
Differences
There are 2 additional lines and one changed line in the new file. After setting up the cell array of defaults, I now check
the varargin cell array for empty arguments. There are ones for which I want to retain the default values. So I find the correct indices
for them, and only overwrite or add the non-empty values from varargin.
Any More Thoughts on Defaults?
If you have any more thoughts on default values, please post them here.
Last week, I wrote a post on the switch statement in MATLAB. One theme in the comments was about using switch (or not!) for setting default values for input arguments that users didn't initialized. I realized that there is a nice
pattern for setting these values that uses some compact, but readable code.
Sometimes I want to write a function that has some required inputs and some optional trailing arguments. If the arguments
are specified by their order, and not by parameter-value pairs, there is a nice way to accomplish this take advantage of varargin. Note that neither of these methods checks the validity of the overridden elements.
Suppose my function requires the first 2 inputs, but there are 3 others that the user can choose to set, or allow them to
take default values. In this scenario, once they choose to set an optional input, they must set all the optional ones that
precede it in the argument list. Here's an example function header.
dbtype somefun21
1 function y = somefun2(a,b,opt1,opt2,opt3)
Here's another way I could write this function header, using varargin.
dbtype somefun2Alt1
1 function y = somefun2Alt(a,b,varargin)
Setting Default Values
To set default values in somefun2, I could use a switch statement or use if-elseif constructs. Here I chose to use a switch statement.
type somefun2
function y = somefun2(a,b,opt1,opt2,opt3)
% Some function with 2 required inputs, 3 optional.
% Check number of inputs.
if nargin > 5
error('myfuns:somefun2:TooManyInputs', ...
'requires at most 3 optional inputs');
end
% Fill in unset optional values.
switch nargin
case 2
opt1 = eps;
opt2 = 17;
opt3 = @magic;
case 3
opt2 = 17;
opt3 = @magic;
case 4
opt3 = @magic;
end
The code is verbose and, in my opinion, not very pretty. It's also error-prone. If I decide to change a default setting
for one value, I have to update each relevant case. What a drag!
Here's a way to set the defaults in one location and then overwrite the ones the user specified.
type somefun2Alt
function y = somefun2Alt(a,b,varargin)
% Some function that requires 2 inputs and has some optional inputs.
% only want 3 optional inputs at most
numvarargs = length(varargin);
if numvarargs > 3
error('myfuns:somefun2Alt:TooManyInputs', ...
'requires at most 3 optional inputs');
end
% set defaults for optional inputs
optargs = {eps 17 @magic};
% now put these defaults into the valuesToUse cell array,
% and overwrite the ones specified in varargin.
optargs(1:numvarargs) = varargin;
% or ...
% [optargs{1:numvarargs}] = varargin{:};
% Place optional args in memorable variable names
[tol, mynum, func] = optargs{:};
First I place the default optional values into a cell array optargs. I then copy the cells from varargin to the correct cells in optargs and I have overridden the defaults. I have only one place where the default values are set, and the code doesn't grow dramatically
in length as I add additional optional inputs. Finally, I spread out the cells from varargin into well-named individual variables, if I want to. Since each element in a cell array is its own MATLAB array, and MATLAB
has copy-on-write behavior, I am not creating extra copies of these optional inputs (see this post for more information).
Other Methods for Dealing with Optional Inputs
There are other methods for dealing with optional inputs. These include specifying parameter-value pairs. With or without
such pairs, you can use inputParser to help set the values you use. You can specify the optional input in a struct, with fields containing the various options, and you can use cell arrays, but then you have to decide how to structure and
specify the contents.
Your Techniques for Specifying Optional Values
Do you have a technique for specifying optional values that works well for you? Please share with us here.
In C, you often terminate case statements in switch constructs with a break statement to prevent execution of the code in the case that follows. Omitting the break statement causes the code in the following case to execute. The behavior of switch statements in C is sometimes described as fall through.
switch in MATLAB
There are several major differences between switch statements in MATLAB and C.
MATLAB does not have fall through behavior in its switch statements and therefore no break statement is typically required for individual case statementss.
MATLAB case statements have another feature that C and many other languages do not have. The expression following case need not be scalar-valued, but can combine several expressions by placing the acceptable expressions in a cell array.
If you want to specify behavior for switch expressions that don't match any of the case expressions, use the otherwise branch. In C, use default instead.
Expressions for case statements can include strings and numeric expressions, not just integer values. Though you have to be careful using non-integral
numbers because of round-off, this flexibility can make the intent of your code more obvious instead of seeing a list of numbers.
More on switch in MATLAB
When we originally were designing switch for MATLAB, we examined much of our own C source code, and we found a huge majority of the case statements terminated in break, hence the no fall-through behavior of MATLAB. This makes the switch statement in MATLAB equivalent to an if-elseif construct. To achieve fall-through behavior in MATLAB, you can specify all the relevant expressions in one case and then conditionally compute values within that section of code.
Do You Use switch Statements?
I am curious to hear about your experiences using switch statements in MATLAB. Post your thoughts here.
For me, the first step in learning is finding out that some feature or technique exists. This often happens in a seminar
or a discussion, sometimes by reading a paper. For me it helps if there's motivation to give me some context. And it's especially
helpful if there's an example solving a problem similar to one I want to solve!
Second Step - Grab Code and Modify
Having decided it's time to learn about some particular algorithm or some feature, what to do next? For some, the first instinct
will likely be to search the web looking for code or references. If I already have the software, I tend to look in the documentation,
and, if possible, find example code close enough to what I hope to do. I then make a copy and modify from there. But that's
just my learning style.
What's Your MATLAB Learning Style?
I'm guessing since you are reading this blog that reading examples is at least partially help to you. There are some other
resources, beyond all the blogs and code on MATLAB Central.
For people new to MATLAB, you can get a great start here.
How do you learn best? Where do you go for your MATLAB information? Post your thoughts here and please include links to any additional valuable references you use for learning MATLAB.